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ABSTRACT: 

A  unified  heuristic  model  of  fluid  turbulence  is  proposed  which  supplements 
the  basic  equations  of  motion  and  continuity  so  as  to  define  determinate 
solutions  under  a  wide  variety  of  conditions.   A  single  set  of  equations 
applies  to  all  cases.   Only  the  boundary  conditions  change  for  each 
application. 

The  Reynolds  stresses  are  related  to  the  mean  flow  heuristically ,  through 
an  eddy  viscosity.   This  also  involves  the  turbulent  energy,  which  is 
analyzed  with  the  aid  of  empirical  expressions  for  dissipation  and  diffusion. 
A  characteristic  local  length  scale  is  introduced,  based  on  a  generalization 
of  von  Karman's  mixing  length. 

The  model  is  illustrated  by  application  to  a  channel  and  a  jet.   Agreement 
with  experiment  is  satisfactory. 

The  report  also  gives  the  detailed  rationale  which  underlies  the  model, 
including  an  original  analysis  of  the  stress/strain  rate  relation.   Suggestions 
are  made  for  further  research. 
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Tart  I 

The  Heuristic  Model  and  Its  Application 

by 

T.  H.  Gawain,  D.Sc. 

and 

J.  V.  Fritcbet 


1.1  INTRODUCTION 

The  Naval  Radiological  Defense  Laboratory,  prior  to  its  disestablish- 
ment in  June  1969,  had  for  some  time  beer  pursuing  studies  in  numerical 
hydrodynamics  in  the  field  of  underwater  explosions.   It  is  veil  known 
that  fluid  turbulence  forms  an  important  aspect  of  this  problem.   In  fact, 
the  lack  of  an  adequate  theory  of  turbulence  has  been  an  important  limita- 
tion on  the  development  of  numerical  hydrodynamics  for  applications  not 
only  to  the  explosion  problem,  but  also  to  fluid  mechanics  problems  of 
every  type.   There  continues  to  be  an  urgent  need  to  overcome  this  obstacle. 

Consequently,  an  agreement  was  reached  between  NRDL  and  the  Naval 
Postgraduate  School  some  time  ago  for  a  joint  attack  on  the  problems  of 
fluid  turbulence.   Two  main  approaches  were  investigated.   The  first  approach 
emphasized  fundamental  understanding  of  the  phenomena  without  undue  concern 
about  developing  computational  methods  for  application  to  practical  problems. 
The  second  approach  emphasizes  the  development  of  a  feasible,  heuristic 
computational  method  at  the  expense,  if  necessary,  of  analytical  rigor. 
This  report  documents  the  excellent  results  that  have  been  achieved  in 
connection  with  this  latter  goal.   Reference  (1),  although  now  largely 
superseded  by  more  recent  work,  documents  some  earlier  effort  along  the 
other  avenue  of  attack.   Further  work  is  continuing  at  the  Postgraduate 
School . 

Unfortunately,  the  Naval  Radiological  Laboratory  was  disestablished 
in  the  midst  of  this  research.   The  abrupt  cessation  of  NPUL  activities 
made  it  virtually  impossible  to  complete  certain  details  of  the 


mathematical  model  in  a  satisfactory  manner.   It  was  therefore  decided  to 
hold  up  distribution  of  a  final  report  until  such  time  as  these  additional 
points  could  be  worked  out.   Meanwhile,  the  press  of  other  sponsored  research 
compelled  the  investigators  to  lay  this  task  aside  for  a  considerable  time. 

Fortunately,  it  was  eventually  possible  to  complete  and  revise  the 
details  mentioned  above.   Consequently,  the  present  report,  which  describes 
the  final  model,  is  now  being  released  for  publication  at  last.   Although 
the  major  part  of  the  work  described  was  accomplished  some  time  ago,  the 
model  itself  should  continue  to  be  of  scientific  interest  for  some  time, 
owing  to  the  basic  and  long  term  nature  of  the  problem  with  which  it  deals. 

As  an  interesting  historic  note,  it  seems  worthwhile  to  mention  that 
the  heuristic  model  developed  in  this  report  resembles  in  some  respects 
the  approach  first  suggested  by  Prandtl  in  19A5,  Ref.  (2).   It  also  has 
some  similarities  to,  and  some  differences  from,  the  more  recent  work  of 
Harlow  and  his  associates,  References  (3)  and  (A). 

In  incompressible  flow,  whether  laminar  or  turbulent,  the  equations 
of  momentum  and  continuity,  along  with  the  boundary  conditions,  theoretically 
suffice  to  establish  completely  the  exact  fluid  motion.   If  the  flow  happens 
to  be  turbulent,  however,  the  actual  detailed  motion  becomes  so  complex 
that,  although  it  is  in  principle  determinate,  its  actual  calculation 
would  impose  an  overwhelming  computational  burden. 

The  usual  way  around  this  difficulty  is  to  average  the  equations  of 
momentum  and  continuity.   This  results,  of  course,  in  an  enormous 
simplification.   Unfortunately,  it  also  involves  a  significant  and  irre- 
trievalbe  loss  of  essential  information.   Consequently,  owing  to  the 
presence  of  the  unknown  Reynolds  stresses  which  are  created  by  this  aver- 
aging process,  the  averaged  equations  of  momentum  and  continuity  do  not 
in  themselves  comprise  a  determinate  set. 


In  order  to  define  a  determinate  solution,  additional  relations  are 
required  to  fix  the  unknown  Reynolds  stresses.   Unfortunately,  these 
supplementary  relations  cannot  he  estahlished  from  the  original  equations 
by  any  purely  deductive  process.   For  this  purpose,  supplementary  empirical 
hypotheses  are  an  unavoidable  necessity. 

From  another  viewpoint,  it  may  be  stated  that  the  averaged  equations 
of  motion  show  the  effect  of  the  Reynolds  stresses  upon  the  mean  flow. 
However,  the  reciprocal  effect  of  the  mean  flow  upon  the  Reynolds  stresses 
is  lost  in  the  averaging  process,   hence  some  adequate  hypothesis  roust  be 
found  for  representing  this  relation,  at  least  approximately. 

For  this  purpose,  a  heuristic  approach  which  seems  plausible  is  to 
postulate  a  relation  between  the  Reynolds  stresses  and  the  mean  flow  which 
is  analogous  to  the  relation  that  is  known  to  govern  the  viscous  stresses. 
The  analogue  of  the  ordinary  molecular  kinematic  viscosity  is  the  so- 
called  eddy  kinematic  viscosity.   The  problem  becomes,  therefore,  to 
determine  empirically  the  general  law  which  governs  this  mean  effective 
eddy  viscosity  at  every  space/time  point  in  the  flow  field. 

The  eddy  viscosity  presumably  depends  on  a  number  of  variables,  one 
of  the  most  important  of  which  is  the  local  kinetic  energy  of  turbulence. 
Therefore,  it  becomes  necessary  to  find  the  space/time  distribution  of  the 
turbulent  energy.   Fortunately,  the  governing  energy  equation  can  be 
deduced  rigorously  from  the  original  equations  of  motion.   However,  the 
energy  equation  itself  introduces  two  additional  unknowns  which  can  only 
be  approximated  in  the  same  heuristic  and  empirical  fashion  as  the  eddy 
viscosity  itself.   The  additional  unknowns  are  the  rate  of  dissipation  of 
turbulent  energy  into  heat,  and  the  rate  of  turbulent  diffusion  of  energy. 

Theory  and  experiment  both  show  that  the  eddy  viscosity,  and  the 
dissipation  and  diffusion  functions  as  well,  depend  not  only  on  the 
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turbulent  energy  itself,  but  also  on  a  local  length  scale  parameter  which 
can  be  associated  with  each  space/time  point  in  the  flow  field.   Von 
Karman  was  perhaps  the  first  to  point  out  how  a  physically  meaningful 
characteristic  length  can  be  defined  in  terms  of  local  space  derivatives 
of  the  mean  velocity  at  any  point  in  the  flow.   In  the  present  paper,  the 
original  approach  of  von  Karman  is  further  developed  and  refined.   It  no« 
takes  into  account  not  only  the  velocity  derivatives  at  the  designated 
point  itself,  but  also  the  values  in  the  general  vicinity  of  the  point. 

By  employing  dimensional  analysis,  and  by  applying  the  available 
experimantal  data,  we  finally  obtain  three  empirical  expressions  which 
determine  to  a  reasonable  approximation  the  eddy  viscosity,  the  heat 
dissipation  and  the  turbulent  diffusion,  respectively.   These  expressions 
also  involve  the  turbulent  energy,  the  local  length  parameter,  and  the 
distance  to  the  nearest  fixed  wall  (if  any) .   Of  course,  these  empirical 
expressions  are  amenable  to  further  investigation  and  development. 

In  this  way  a  single  consistent  and  determinate  set  of  equations  is 
established  which  applies  in  principle  to  any  incompressible  turbulent 
flow  field.   Only  the  boundary  conditions  differ  for  each  specific 
application . 

The  application  of  this  heuristic  model  is  shown  for  two  widely 
different  cases,  namely,  a  two-dimensional  duct,  and  an  axi-symmetric  jet. 
The  distributions  of  energy  and  Peynclds  stress  are  computed  for  each 
case  and  compared  with  experimental  measurements.   The  agreement  is 
satisfactory. 

These  results,  although  limited  in  extent,  strongly  suggest  that  the 
proposed  heuristic  model  of  fluid  turbulence  is  basically  correct,  at 
least  in  its  main  essentials. 


This  report  is  divided  into  two  parts.   Part  I  summarizes  the  final 
heuristic  model  and  illustrates  its  application  to  two  typical  cases, 
namely,  a  channel  flow  and  a  free  jet.   Part  II  fives  the  corresponding 
detailed  background,  theory  and  rationale.   It  tells  vhy  the  final  model 
is  just  what  it  is.   It  shows  bow  the  various  parameters  of  the  model  are 
related  to  fundamental  aspects  of  turbulent  correlation  and  spectra.   An 
original  analysis  is  presented  of  the  general  tensor  relation  between 
stress  and  strain  rate.   This  provides  the  basis  for  judging  both  the 
usefulness  and  the  limitations  of  the  eddy  viscosity  concept  as  employed 
in  the  heuristic  model.   The  energy  equation  is  derived  from  first 
principles,  then  transformed  so  as  to  be  practically  useful.   This  part  of 
the  report  also  raises  various  questions  of  a  fundamental  nature  which 
could  be  fruitful  for  future  research. 


1.2  SYMBOLS 

Fef.  Eos. 
a  Local  dimensionless  length  scale  in  jet     1.6-31 

a1  Overall  dimensionless  length  scale  of       1.6-28 

cross-section  in  jet 

A  Cross-sectional  area 

b  Channel  half  width, or  jet  reference  length   1.6-1 


AlA2A3| 

B1B2B3^ 


ij 


Functions  in  momentum  equations  for  jet     1.6-9  thru  -14 


2    2    2 
Exponents  in  definitions  of  I  ,  J  ,  X  ,     2.6-14  thru  -19 

and  t 

Arbitrary  constants  in  definitions  of       2.6-8,-9 
A  2  and  i 

D  Determinant.   Also  diameter  of  jet  2.6-16 

E  Mean  kinetic  energy  of  turbulence  per       1.3-16,  2.3-5 

unit  mass 

E'  Instantaneous  kinetic  energy  of  turbulence  A-12 

per  unit  mass 

E  Rate  of  dissipation  of  turbulent  energy     2.8-8 

per  unit  mass  into  heat. 

f  Gaussian  function  1.6-6 

f  Ratio  of  local  length  scale  of  mean  flow    2.5-1 

to  correlation  length 

f  Ratio  of  local  time  scale  of  mean  flow      2.5-2 

to  correlation  time 

f..  Non-isotropic  tensor  in  generalized  stress/  2.7-16 


strain  rate  law 


f *  .  Denotes  values  of  f..  for  principal  axes 

n  11 

J  of  strain  rate 


Dimensionless  stream  function  of  jet        1.6-1 

fdF 
Denotes  derivative  — 

[in 

Functions  in  energy  equation  for  two-       1.5-25,  -26-27 
dimensional  channel 

Parameter  in  stress/strain  rate  law  2.7-28 

Functions  in  analysis  of  jet  1.6-25,  -26,  -27 

Integral  of  weighting  function  1.3-11 

Integral  of  spatial  weighting  function  2.5-9 

Integral  of  temporal  weighting  function  2.5-10 

Tensor  indices 

Characteristic  integrals  for  a  point  1.3-13  and  -14 

Von  Kantian' s  mixing   length  1.5-7,    -8 

Superscript  indicating  realization  number 

Reynolds  number 

Kinematic  pressure  (pressure  divided  by 
density)  of  turbulent  fluctuations 

P  Kinematic  pressure  (pressure  divided  by 

density)  of  mean  flow 

q  Root  mean  square  turbulent  velocity         1.3-16 

fluctuation 

r,  z,  G      Cylindrical  coordinates:   radial,  axial, 
and  angular,  respectively 

R .  .  Correlation  tensor  2.3-6 

ij 

(S) 
K:/  Spatial  correlation  tensor  2.3-9 

ij 

(T) 
R^  Temporal  correlation  tensor  2.3-9 

s  Exponential  constant  in  Gaussian  velocity   1.6-6 

distribution  for  jet  (s  =  0.102) 

t  Time 


F 

F' 

F,    G,   H. 

g* 

Gr  G2' 

G3 

H 

„«> 

H(T) 

i»    j,   k, 

1 

2        2 
I    ,   J 

I 

(n) 

NR 

P 

At  Time  increment                            2.3-1 

T  Semi-interval  for  time  averaging           2.2-4 

t..  Viscous  stress  tensor  for  turbulent         A-16 

~*  fluctuations 

T..  Viscous  stress  tensor  for  mean  flow 

T'  Viscous  stress  tensor  for  mean  flov  plus    A-2 
i  i 

J  turbulent  fluctuation 

u.  Turbulent  velocity  fluctuation             2.2-2 
i 

U.  Mean  velocity                            2.2-1 
l 

U!  Velocity  of  mean  flow  plus  turbulent 
i 

fluctuation 

U,  V  Radial  and  axial  components,  respectively, 
of  axi-symmetric  mean  flov 

U  Reference  velocity  of  turbulent  jet         1.6-1 

o  J 

-u  .u  .^ 

-1  \                Kinematic  Reynolds  stress  tensor  (stress    2.7-1 

-uv  J                 divided  by  density) 

v*  Friction  velocity                         1.5-2 

dv'  Element  of  volume 

W  Weighting  function                         1.3-10 

(S) 

W  Spatial  weighting  function                 2.5-4 

(T) 

W  Temporal  weighting  function                2.5.5 

w  Normalized  weighting  function              1.3-12 

(S) 

w  Normalized  spatial  weighting  function       2.5-12 


t  (T) 

w 


—  i 

X 


Normalized  temporal  weighting  function      2.5-13 

Position  vector  of  arbitrary  fixed 
reference  point 


Position  vector  of  variable  point 
Ax  Displacement  vector  1.3-9 


*i     1 

>    Cartesian  coordinates 
x,  y,  zj 

y  Distance  to  nearest  wall ,  or  distance 

from  centerline  of  jet  or  center  plane 
of  channel 

z  Axial  coordinate.   Also  nondimensional 

wall  distance  exponent  1.5-23 

a  Dimensionless  eddy  viscosity  coefficient    1.3-5,-18 

P,  $'        Dimensionless  heat  dissipation  coefficients  1.3-19,2.8-10 

Y  Dimensionless  diffusion  coefficient         1.3-20,2.8-21 

Y. .  Strain  rate  tensor  of  turbulent  A-17 

J  fluctuations 

P.,  Strain  rate  tensor  of  mean  flow  1.3-4 

ij 

T*  ,T*  ,Ti^  Strain  rate  tensor  components  for  principal 
axes  of  strain  rate 

6..  The  Kronecker  delta  tensor,  with  value  unity 

J  for  equal  indices,  but  zero  otherwise 

e  Effective  eddy  kinematic  viscosity         2.7-22 

e  Eddy  viscosity  in  centerline  of  jet         1.6-16 

e,,  .£..,,    Tensor  versions  of  eddy  viscosity  2.7-13,-14 

jk'  ljkl 

e"  "Pressure  diffusion'  coefficient  2.8-17 

e'  Kinetic  energy  diffusion  coefficient        2,8-18 

z,  Radial  coordinate  B- 11 

n  Dimensionless  distance  from  center  plane    1.5-1 

of  two-dimensional  channel,  or  from 
centerline  of  jet 

0,  61        Angular  coordinate 

6*  Strain  rate  parameter  in  stress/strain  rate  2.7-26 

law 

k  Constant  in  von  Karman's  mixing  length      1.5-8 

theory 
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k,  k.  Wave  number  vector  1.5-8 

A  Local  length  scale  of  mean  flow  1.3-15 

A  Overall  length  scale  1.3-8 

A*  Correlation  distance  2.3-11 

A  Dissipation  length  2.8-8 

v  Molecular  kinematic  viscosity 

§  Separation  distance  in  meridional  plane  B-6 

I  Summation  sign 

t  Local  time  scale  of  mean  flow  2.6-18 

t*  Correlation  time  2.3-12 

i..  Kinematic  Reynolds  stress  tensor  (stress    2.7-1 
i  i 

divided  by  density) 


ij 


,-2 


Distortional  kinematic  stress  tensor        2.7-4 


t..  Residual  non-isotropic  kinematic  stress     2.7-19 

J  tensor 

(t!. )        Viscous  distortional  stress  tensor  2.7-10 


Invariant  of  the  tensor  t'  2.7-5 

2 
t"  Invariant  of  the  tensor  t  , ,  2.7-20 

<J> .  .  Fourier  transform  of  the  correlation        2.A-2 

J  tensor  R . . 

ty  Stream  function  for  jet.   Also  auxiliary    1.6-1 

weighting  function  for  axi-symmetric  flows  B-  13 

i>*  Stress  parameter  in  stress/strain  rate  law  2.7-29 

w  Angular  frequency  or  temporal  wave  number 

ft  Invariant  of  strain  rate  tensor  1.3-6 

representing  generalized  shearing  rate 

ft1  Absolute  magnitude  of  the  gradient  of  ft     1.3-7,2.6-5 

Note:   In  the  jet  analysis,  subscript  1  on  any  quantity   1.6-33 
denotes  the  values  at  z  =  1.   It  is  therefore 
always  a  function  of  n  only. 
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1.3  SlfMMAFY  OF  BASIC  EQUATIONS  OF  THF  FFUP.ISTIC  MODEL 


The  purpose  of  this  section  is  to  summarize  the  basic  eouations 
which, taken  together,  constitute  the  proposed  heuristic  model  of  fluid 
turbulence.   A  few  explanatory  comments  arc  offered  to  identify  these 
various  equations  and  to  show  their  general  relation  to  the  overall  theory, 
but  no  detailed  justification  is  given  in  this  section;  such  details  are 
included  in  other  parts  of  this  report.   Also,  the  eouations  in  this 
section  are  stated  in  three  dimensional  cartesian  form.   The  reductions 
to  forms  specifically  applicable  to  plane  or  axi-symmetrical  mean  flows 
are  treated  elsewhere. 

An  incompressible  mean  flow  field  must  everywhere  satisfy  the 


continuity  condition 

3U, 

=  C 


3x 


(1.3-1) 


and  the  ensemble  averaged  Navier  Stokes  equations  of  motion 


3U. 

l 


3t 


U.U. 
3x. 1  l  jj 


=  v 


2 

31) 


3x .  3x  . 
J   J 


3P 


i33V 


>\? 


3x 


J 


(-2) 


as  well  as  the  particular  boundary  conditions  appropriate  to  the  specific 
flow  field  under  consideration. 

However,  because  of  the  occurence  of  the  unknown  Feynolds  stress 
x   ,  Fqs.  (-1)  and  (-2)  do  not  in  themselves  constitute  a  determinate  set 
of  equations.   In  order  to  define  a  determinate  solution,  additional  relations 
are  required.   Purely  deductive  theory  alone  does  not  provide  these  relations. 
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Hence  appropriate  supplementary  hypotheses  must  be  tentatively  postulated, 
subject  to  eventual  verification,  modification  or  disproof  by  experimental 
evidence. 

For  the  reason  discussed  at  length  elsewhere,  we  postulate  that 
the  Reynolds  stresses  can  be  adeauately  related  to  the  strain  rates  of  the 
mean  flow  through  a  law  which  has  been  deliberately  somewhat  oversimplified 
to  the  form 

t..  =  -  i  qV  .  +  el\  .  (-3) 

2 
In  this  equation,  q   is  the  mean  square  velocity  of  the  turbulent 

fluctuations,  e  is  the  so-called  eddy  kinematic  viscosity,  and  r  .  is  the 

strain  rate  of  the  mean  flow  as  defined  by  the  expression 

(-4) 


r. . 


3U .    3U . 


3x .    3x . 


Unfortunately,  the  eddy  viscosity  e  as  defined  above  is  not  a  simple 
property  of  the  fluid  itself,  but  rather  a  complex  property  of  the  turbulent 
flow  field.   Physical  and  dimensional  reasoning  leads  to  a  supplemantary 
postulate,  namely  that  c    is  expressible  in  the  following  form. 


_€ 


=  a  (-5) 


Here  A  represents  a  local  length  scale,  defined  below,  which  varies 
from  point  to  point  in  the  field,  but  which  is  everywhere  related  in  a 
consistent  way  to  the  local  mean  velocity  distribution.   The  dimensionless 
coefficient  a  is  a  slowly  varying  universal  function  not  predictable  from 
theory,  but  determinable  from  experimental  data.   This  function  is  specified 
later  in  Eq.  (-18) 
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To  establish  the  quantity  A  at  any  point  x,  we  proceed  as  follows. 
We  first  compute  the  following  mean  flow  variables  throughout  the  field, 
namely, 

a2  -  -i  r  r 

2  'ij'ij 


(-6) 


2    1 


9fi 


3x 


an' 


3x 


(-7) 


Then  we  compute  the  reference  quantity 


^2  = 


/fi4dv 


/(^')2dv 


(-8) 


where  the  volume  integrals  extend  over  the  entire  field.  Of  course  A   is 
therefore  a  characteristic  constant  for  the  field  as  a  whole  and  is 
independent  of  x,  but  it  may  vary  as  a  function  of  time  if  the  mean  flow 
be  unsteady. 

The  calculation  continues  as  follows.   All  integrals  extend  over 
the  entire  field  of  flow.   Let 

x  =  coordinate  of  the  fixed  point  at  which  A  is  being  evaluated 
5c'  =  coordinate  of  variable  point  at  which  integrands  below  are 

evaluated 
Ax  =  ( x'  -  x)  -  separation  variable  (-9) 


Ax-Ax 

"  ",  2 


W  =   e 


=  weighting   function 


(-10) 


H  =  /Wdv'     =  normalizing  integral 


w  = 


=  normalized  weighting  function 


(-11) 
(-32) 
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2      A   , 
T   =  fwQ.   dv  =  first  characteristic  integral  (-13) 

2  o 

J   =  /w(f2fi')  dv  =  second  characteristic  integral  (-14) 

2      2 
Of  course  the  integrals  T   and  J   are  functions  of  x  and,  if  the 

mean  flow  is  unsteady,  of  time  t  as  well. 

The  required  length  scale  A  at  any  point  x  is  finally  given  by  the 

ratio 

2    I2 
AZ  =  ~  (-15) 


The  above  calculation  sequence  establishes  not  only  the  desired 
local  length  scale  A  but  also  the  useful  auxiliary  quantity  J.   It  turns 
out  that  J  is  needed  below  in  connection  with  the  energy  eciuation.   The 
rationale  involved  is  fully  explained  elsewhere  in  this  report. 

If  the  mean  flow  happens  to  be  either  plane  or  axi-symmetri c ,  all 
of  the  above  volume  integrals  can  be  reduced  to  corresponding  surface 
integrals.   The  reduction  for  the  axi-symmetric  case  is  rather  interesting; 
it  is  given  in  Appendix  B. 

If  the  flow  field  happens  to  be  self-similar,  the  above  calculations 
need  be  made  on  only  one  representative  cross-section.   Corresponding 
results  for  other  cross-sections  are  then  found  simply  by  applying  the 
proper  scaling  rules.   If  the  flow  field  is  not  self-similar,  all  points 
must  be  computed  independently.   If  the  flow  field  is  steady,  the  calcula- 
tions can  be  made  once  and  for  all.   If  it  is  unsteady,  they  must,  of 
course,  be  repeated  at  appropriate  intervals  of  time. 

Once  the  distribution  of  A  has  finally  been  established,  it  is  then 
necessary  to  determine  the  distribution  of  the  turbulent  energy  E,  or 
equivalently ,  of  the  root  mean  souare  turbulent  velocity 
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q  =y2E 
This  is  governed  by  the  energy  equation 


(-16) 


=  en2   -   B(2E)7/6J1/3  -  ^-[D.E]  +  -^- 

Sx.   i     3x. 
1  1 


Y£ 


3E 


L3x1 


J 


(-17) 


Of  course  the  eddy  viscosity  e  in  the  energy  equation  is  the  sane 
as  that  earlier  defined  in  Eq .  (-5)  in  connection  with  the  stress/strain 
rate  law.   If  desired,  Eq .  (-5)  may  be  used  to  eliminate  e  from  the 
energy  equation.   In  any  case,  the  coefficient  a  which  governs  e  is  fixed 
by  Eq.  (-18)  below. 

Each  term  in  Eq .  (-17)  has  dimensions  of  time  rate  of  change  of 
energy  per  unit  mass.   The  four  terms  on  the  right  represent,  respectively, 
generation  of  turbulent  energy,  dissipation  to  heat,  convection  by  the 
mean  flow,  and  turbulent  diffusion.   Molecular  diffusion  is  negligible. 

The  dimensionless  coefficients  3  and  y  which  occur  in  the  energy 
equation,  like  the  dimensionless  coefficient  a  introduced  earlier,  are 
assumed  to  be  well  behaved  universal  functions,   however,  they  are  not 
determinable  from  theory,  but  must  be  established  empirically. 

On  the  basis  of  the  work  done  so  far,  the  following  expressions 
seem  to  provide  satisfactory  agreement  with  the  available  data. 


a  =  0.065<1  +  e 


-•?-» 


(-18) 


4  -  3 . 7  a  +  e 


.I*.: 


Y  =  1.4  -  O.Ae 


-'f- 


(-19) 


(-20) 
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In  these  expressions,  y  is  the  distance  to  the  nearest  fixed 
boundary.   As  the  boundary  is  approached  y  approaches  zero,  but  it  turns 
out  that  under  these  conditions,  the  length  parameter  A  also  approaches 

y 

zero  simultaneously,  so  that  the  ratio  *■  remains  finite.   In  fact,  it 


happens  that  at  the  wall  itself,  -*•- 


equals  unity  (not  zero).   In  any  case, 
it  is  stipulated  in  connection  vith  the  above  forrrulas,  that  \~-\    shall 
arbitrarily  be  assigned  a  lower  limit  of  unity.   On  the  other  hand,  far 
from  any  fixed  boundary,  or  in  the  absence  of  fixed  boundaries,  ve  set 
y  =  <»,  whereupon  a,  f,  and  y  reduce  to  three  simple  constants. 

The  proposed  heuristic  model  of  fluid  turbulence  consists,  therefore, 
of  the  foregoing  equations  (-1)  through  (-20)  inclusive,  plus  the  particular 
boundary  conditions  appropriate  to  the  specific  problem  under  consideration. 
This  system  of  eauations  theoretically  defines  a  determinate  solution  for 
the  general  case  of  inhomogeneous  and  non-stationary  fluid  turbulence.   It 
fixes  the  mean  flow  field,  the  distribution  of  turbulent  energy  and  the 
Reynolds  stresses. 

The.  results  obtained  so  far  tend  to  support  the  proposed  model,  but 
further  work  is  necessary  to  ascertain  more  fully  its  adenuacy  and  range  of 
applicability.   Various  features  of  the  model  arc  readily  amenable  to 
modification  if  necessary,  so  as  to  improve  its  ranpe  and  accuracy. 

An  important  merit  of  the  model,  if  it  continues  to  prove  successful, 
is  that  it  reduces  a  wide  range  of  dissimilar  flow  geometries  and  flox* 
conditions  to  a  single  unified,  theory. 

1.4  CEKEFAI  METHOD  OF  SOLUTION 

Although  the  equations  (-1)  through  (-20)  of  the  previous  section, 

along  with  the  appropriate,  boundary  conditions,  theoretically  suffice  to 
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determine  all  characteristics  of  the  flow  field,  the  actual  solution  of 
this  system  of  equations  for  any  non-trivial  case  is  by  no  means  simple. 
Let  us  consider  first  the  general  case  of  an  unsteady  mean  flov,  with  the 
initial  velocity  distribution  of  time  t  =  0  arbitrarily  prescribed. 

Of  course,  the  prescribed  initial  velocity  distribution  must  satisfy 
the  continuity  condition,  Eq.  (-1).   This  offers  no  particular  difficulty. 
In  the  usual  case  in  which  the  mean  flov  is  either  plane  or  axi-symmetri c, 
the  continuity  condition  may  be  satisfied  identically  by  expressing  the 
velocity  components  in  terms  of  a  stream  function.   Let  us  suppose  that 
this  has  been  done,  and  that  the  corresponding  mean  flow  velocity  components 
at  time  t=0  have  been  calculated.   The  surseauent  calculation  procedure  may 
now  be  described  in  somewhat  general  terms  as  follows. 

The  next  step  is  to  perform  the  sequence  of  calculations  summarized 
in  Eqs.  (-4)  through  (-15).   In  this  way  is  obtained  the  local  length 
scale  >  at  every  point  in  the  field,  and  the  ouantity  J. 

Next  we  must  solve  the  energy  equation  (-17)  for  the  mean  souare 

turbulent  velocity  o  at  every  field  point.   For  this  purpose  ve  substitute 

Eq.  (-5)  into  (-17)  and  use  Eqs.  (-18),  (-19),  and  (-20)  to  evaluate  the 

three  empirical  coefficients  a,  £,  and  y.      (In  some  cases  the  actual 

2 
distribution  of  the  turbulent  energy  E  =  -^  may  be  known  from  experimental 

measurement.   If  so,  we  have  at  this  point  one  important  comparison  between 

theory  and  experiment.) 

Once  estimates  of  A  and  n  are  available  from  the  foregoing  calculations, 

it  becomes  a  simple  matter  to  solve  Eq.  (-5)  for  the  eddy  viscosity  e,  and 

Eq .  (-3)  for  the  Reynolds  stresses. 
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The  Reynolds  stresses  obtained  in  this  way  may  now  be  substituted 
back  into  the  Navier  Stokes  relations,  Fes.  (-2).   By  taking  the  divergence 
of  Cqs.  (-2),  we  find  that 


3  P 


[l.b.  +  t.  . 


3x . 3x  .   i  j    i  i 
i  J 


3x . 3x  . 

1  *  r 

Since  mean  velocity  and  stress  are  at  this  point  known,  the  al ove  result 

can  be  solved  for  the  pressure  distribution,  subject  to  the  specified 

boundary  conditions  on  pressure.   This  part  of  the  solution  involves  a 

lengthy  iteration  procedure,  but  is  straightf orward  in  princinle. 

Once  l . ,  x...  and  P  are  all  established,  we  can  revert  to  En.  (-2) 

/-ST  ,i 
and  solve  for  the  time  rate  of  chanpe  of  velocity 


(1.4-1) 


[3t  j 


at  every  field 


point.   By  integrating,  forward  through  time  over  a  small  interval,  we 
can  find  the  new  mean  velocity  distrit  ution  at  the  later  time. 

The  above  procedure  may  be  repeated  as  often  as  desired,  thus 
tracing  the  gradual  evolution  of  the  flow  field  over  time.   If  the  boundary 
conditions  happen  to  be  stationary  in  time,  the  entire  flow  field  wil 1 
evolve  toward  a  stationary  state,  thus  defining  the  solution  of  the 
corresponding  steady  state  problem. 

The  above  process  of  reaching  a  steady  state  asymptotically  in 
time  can  be  quite  expensive  in  terms  of  computational  effort  and  time 
required.   It  is  often  possible  to  shorten  or  abridge  the  above  calculation 
in  various  ways,  depending  on  the  situation.   For  example,  instead  of 
substituting  the  Reynolds  stresses  predicted  by  F.o .  (-3)  into  the  Navier 
Stokes  lav;  Eq.  (-2),  we  can  solve  the  latter  for  the  Reynolds  stresses 
'required"  in  order  to  maintain  the  prescribed  mean  motion  in  a  steady 
state.   These  "required"  stresses  can  then  be  compared  with  the  stresses 
'predicted"  from  Eq .  (-3).   If  the  agreement  between  the  two  sets  of 
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stresses  is  deemed  satisfactory,  this  can  he  taken  as  a  verification  of  the 
adequacy  cf  the  given  velocity  distribution  as  a  estimate  of  the  steady 
state.   This  is  the  procedure  that  has  been  followed  in  our  sample  problems, 
and  it  has  proved  entirely  satisfactory  for  the  intended  purpose.   In 
this  case  the  purpose  has  been  to  demonstrate  the  degree  of  agreement  of 
our  heuristic  model  with  experimental  observations.   The  agreement  is 
considered  satisfactory  for  both  cases  analyzed. 

1.5  APPLICATION  TO  A  TWO-DIMENSIONAL  CHANNEL 

Consider  the  steady,  incompressible,  turbulent  flow  through  a  two- 
dimensional  channel  of  uniform  height  2b  as  shown  in  Fig.  1.5-1.   Let  y 
be  distance  from  the  channel  center-line.   Choosing  b  as  the  reference  unit 
of  length,  the  corresponding  dimensionless  distance  becomes 

I  =  n  (1.5-1) 

Let  the  kinematic  shear  stress  at  the  wall  be 

x1  =  v*2  (-2) 

where  v*   is  the  so-called  friction  velocity.   Without  loss  of  generality, 
we  can  always  choose  reference  units  of  time  such  that 

t1  =  v*2  =  1  (-3) 

The  kinematic  shear  stress  t  varies  linearly  across  the  channel  so 
that  at  any  location  n  we  have 

x  =  n  (-4) 

According  to  our  present  notation,  if  I'  is  the  mean  velocity  at  any 
n ,  then 
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FIGURE  1.5-1    FLOW  IN  A  TWO-DIMENSIONAL  CHANNEL 
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dn 


-  U'  -  fi 


dV 


dn 


=  u,:  =  n' 


(-5) 


(-6) 


According  to  the  von  Karman  mixing  length  theory 
2  2 

t  =  i  n 

where  the  mixing  length  I   is  given,  in  our  notation,  by 
„2    2  ft2 


(-7) 


n 


,2 


(-8) 


and  where  k  is  a  constant  which  remains  to  be  determined. 

Upon  combining  Eqs.  (-4),  (-7),  and  (-8),  taking  the  square  root 
of  both  sides  ard  rearranging,  we  obtain  the  basic  relation 


.2 


(-9) 


This  equation  is  readily  integrable  and  gives  the  well  known  von 
Karman  velocity  profile  which  in  the  present  notation  becomes 


IL  =  i  =  u  +  -  {£n(l  -  /n)  +  /n"} 

V*  OK 


(-10) 

It  is  known  that  this  results  gives  satisfactory  agreement  with  experiment 
if  we  place 

k  =  0.36  (-11) 

The  resulting  quantities  of  interest  in  the  present  analysis  are 


n  -  V 


2k(1  -  /n") 


and 


C  -  U"  «  - 


Ak  n(l  -/n) 


(-12) 


(-13) 
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We  next  require  the  'averaging  distance1'  A   defined  by  the  equation 


f 


^2     .  .    '0 

A,  =  £im 


Vdn 


1     ,  f1 
V1  J  o 


(-14) 


(QQ')2dn 


Upon  substituting  expressions  (-12)  and  (-13)  into  (-14) ,  we  find 
that  in  the  limit  as  n,  approaches  unity,  both  integrals  approach  infinity, 
but  the  above  ratio  of  these  integrals  approaches  zero.   Hence  we  obtain 
the  unusual  result  that,  in  this  case, 


Now  the  local  length  scale  A  may  be  computed  from  the  relation 
■1   Ax*Ax 


e   A 


A   =  I  im 
A  -K) 


1  nW 


1  A  x«  A  x 


e   A 


1  (fifi')2dn' 


(-15) 


(-16) 


This  result  is  remarkable!   It  shows  that  our  characteristic 
length  A  is  in  this  case  proportional  to  the  classical  von  Karman  mixing 
length  I.    Form  Eqs .  (-8),  (-12),  (-13)  and  (-16)  we  find  that 


Similarly,  because  of  the  exceptional  result  (-15)  ,  the  integral 


for  J   reduces  simply  to 

(1 


A  x«A  x 


e    l      (M')2dn' 


J   =   £im 
A  +0 


o 


1  _  Ax'Ax 
? 

Xl 

e       dn' 


=  (f^')2 


(-17) 


(-18) 
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Therefore  again  using  (-12)  and  (-13)  we  find  that 


j  =  fin'  = 


8K2/n-(l  -  /n)3 


(-19) 


The  three  empirical  functions  a,  $,  and  y   may  now  be  expressed  in 

the  form 

2 
=  a  =  0.065{1  +  e  Z  } 


A  ^21. 


(-20) 


and 


where 


=  3.7(1  +  e  Z  } 


=  1.4  -  0.4  e 


-z 


2    if  1    .12 

Z  =  aIt^  "  \ 

With  the  functions  ft,  A,  J,  a,  $,  and  y  now  known,  we  can  proceed 
to  solve  the  energy  equation.  The  unsteady  term  and  the  convective  term 
both  vanish  from  this  equation  in  the  present  instance  and  we  obtain 


(-21) 


(-22) 


(-23) 


a\J2YS)2  -  3(2E)7/6J1/3  + -Ir- 
an 


^ 


=  0 


(-24) 


Energy  E  is  now  the  only  unknown.   To  simplify  the  equation  further, 
we  introduce  the  following  notation. 


a\tl       =  F(n) 

1/3 
BJ  /J  =  G(n) 


(-25) 


(-26) 

<xy\       =  H(n)  (-27) 

Of  course,  these  three  functions  are  now  known.   Therefore  the  energy 
equation  reduces  to 
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dn 


dn 


-  G(2E)?/6  +  F/2E  =  0 


(-28) 


This  equation  cannot  be  solved  analytically,  but  it  can  easily  be 
expanded  and  rearranged  for  numerical  solution  by  finite  differences. 
The  desired  boundary  conditions  are 

At   n  =  0    (4^1   =  E'  =  0 


(-29) 


At   n  =  1   E.  -  0 


(-30) 
This  last  condition  cannot  be  satisfied  directly.   However,  we  can  integrate 


from  the  centerline  out,  using  a  trial  value  of  the  centerline  energy  E  . 

By  trial  and  error,  a  value  of  F   can  then  be  found  which  will  in  fact 

o 

yield  the  required  condition  E..  =  0  at  the  wall. 

This  is  the  method  which  was  used  to  compute  the  energy  distribution 
shown  by  the  solid  line,  in  Fig.  1.5-2.   It  may  be  seen  that  the  curve 
obtained  in  this  way  is  in  reasonably  good  agreement  with  the  experimental 
data  of  Feichart  and  of  Laufer,  Refs.  5  and  6. 

Once  the  energy  distribution  is  known,  the  dimensionless  Feynolds 
stresses  may  be  computed  from  the  relation 

-uv  =  t  =  aX,/2Efi  (-31) 

Except  very  near  the  walls,  the  above  Feynolds  stress  should  agree 

with  the  total  shear  stress  defined  by  Eq .  (-4). 

If  we  wish  to  switch  back  to  absolute  units,  Eqs.  (-4)  and  (-31) 

reduce  respectively  to 

(-32) 


1 =  1 

.,2   b 


Total  Stress 


uv  A 

~2=  a\b 


fv/2E 


bfi 


lv~) 


Reynolds   Stress 


(-33) 
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In  Fig.  1.5-3,  the  solid  line  shows  the  above  heuristic  estimate  of 
the  Reynolds  stress,  computed  for  a  Reynolds  number  of  10,000.   The  dashed 
straight  line  represents  the  true  total  stress,  which  should  agree  closely 
with  the  Reynolds  stress  everywhere  except  in  the  immediate  vicinity  of 
the  wall.   The  agreement  is  seen  to  be  good;  the  computed  results  are  on 
the  whole  more  accurate  than  are  the  experimental  measurements  from  Reference 
(6),  as  shown  on  this  figure. 

We  conclude  that  this  example  tends  to  substantiate  the  proposed 
heuristic  model  of  turbulence. 

1.6  APPLICATION  TO  AN  AXI -SYMMETRIC  JET 

Consider  a  free  turbulent  jet  discharging  into  a  quiescent   atmosphere 
as  shown  in  Fig.  1.6-1.   The  radial  and  axial  coordinates  are  r  and  z  and 
the  corresponding  velocity  components  are  U  and  V,  respectively. 

Experimental  results  indicate  that  the  velocity  profiles  of  the 

mean  flox>?  at  various  cross-sections  are  self-similar.   Specifically,  if 

z  be  measured  from  a  suitable  virtual  origin,  the  stream  function  \b   can 

be  reduced  to  the  form 

4>  =   U  bzF(n)  (1.6-1) 

o 

where 

n  =  -  (-2) 

z 

and  where  U  and  b  are  constants  of  the  jet. 
o  J 

Consequently  the  velocity  components  of  the  mean  flow  become 
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FIGURE  1.6-1     THE  AXI-SYMMETRIC    JET 
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where 
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8F 
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The  known  experimental  results  also  Indicate  that  the  generalized 
velocity  distribution  through  the  jet  is  very  well  approximated  by  the 
simple  expression 


—  =  f  =  e 
n 


2 


(-6) 


Where  s  is  a  characteristic  constant  for  a  turbulent  jet.  It  has  the  known 
value  s  =  0.102. 

Furthermore,  it  is  useful  and  convenient  to  adopt  units  of  length 
and  time  such  that,  when  expressed  in  these  units,  U  ■  1  and  b  =  1. 
Of  course,  there  is  no  loss  in  generality  involved  in  this.   It  amounts 
merely  to  non-dimensionalizing  all  quantities  with  respect  to  U  and  b  as 
reference  parameters.   Hence  symbols  like  e,  p,  u,  v  and  so  on  will  now 
represent  the  corresponding  dimensionless  versions  of  kinematic  viscosity, 
kinematic  pressure,  velocity  components  and  the  like. 

The  next  step  in  the  analysis  is  to  combine  the  momentum  eouations 
(-2)  and  the  stress  lav  (-3)  of  section  1.3.   The  results  may  be  transformed 
to  the  z,  n  coordinate  system  defined  above.   The  foregoing  dimensionless 
functions  e  and  p  may  be  introduced  therein,  and  the  results  reduced. 
This  reduction,  although  very  lengthy,  is  straigbtf orward .   In  this  way 
the  following  results  are  obtained. 

np'  +  2p+  A1(n)e'  +  A2(n)e  =  A3<n)  (-7) 

-p'  +  0  +  B1(n)e'  +  B2(n)e  =  B3(n)  (-8) 
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Ihe  A's  and  B's  are  known  functions  defined  as  follows,  where  f  is  the 


Gaussian  function  defined  by  Eq.  (-6). 
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Th.eoretically ,  it  is  possible  to  eliminate  p  and  p'  between  the  two 
momentum  equations  (-7)  and  (-8).   This  yields  a  single  vorticity  equation 
of  second  order  in  c  which  can  be  solved  for  the  unknown  eddy  viscosity 
distribution  c(n). 

However,  it  is  known  from  experiment  that  the  terms  in  p  and  p'  are 
negligible  in  Eqs.  (-7)  and  (-8).   Moreover,  the  terms  in  bq.  (-8)  are 
negligible  in  comparison  with  those  in  Eq .  (-7).   We  can  tale  advantage  of 
this  information  to  simplify  the  solution  for  c(n).   We  simply  discard 
Eq .  (-8),  and  drop  the  terms  in  p  and  p'  from  bq .  (-7).   The  resulting 
equation  is  now  only  of  first  order  in  t.   It  may  be  written  in  the  form 


vdn 
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s2J 
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(-15) 
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Solution  of  this  equation  analytically  is  not  feasible,  but  solution 
by  numerical  integration  is  simple. 

>;ote  first  that  at  the  centerline  n  =  0 ,  this  equation  gives  at  once 


-  1 


=  0.00260 


(-10 


This  result  provides  the  boundary  value  required  for  the  numerical  integration. 

It  is  also  of  interest  to  note  that  for  very  large  values  of  n, 
Eq.  (-15)  simplifies  to  the  following  approximation,  namely, 

2 
•2q" 


dc 

dn 


Id 
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£  =  - 


S 


(-17) 


The  solution  is 


5  .«  2      _    2 

s  +sn  e°°  -5  -sq  , 

e  =  --_—  e     /  q  e    dn 

2  ;  q 


(-18) 


The  true  numerical  solution  approaches   the  above  analytical  curve  asymptoti- 
cally at  large  values  of  q. 

V/ith  the  dimensionless  eddy  viscosity  distribution  c(q)  Vnown,  it 
becomes  a  simple  matter  to  find  the  corresponding  dimensionless  Reynolds 

shear  stress  t    from  the  relation 
rz 


rz 


-UV  =  L 


+ 

(         1 

(-19) 


By  making  use  of  Eqs.  (-3)  through  (-6)  this  result  can  be  reduced  to  the 
following  form  for  computation 


t   =  -uv  = 
rz 


e(n) 


3  +  --- 


q  +  -p       f 


(-20) 


The  shear  stress  computed  from  this  enuation  is  shown  by  the  dashed  line 
in  Fig  (1.6-2). 
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To  form  the  mean  flow  parameters  required  for  the  heuristic  model, 
we  note  that,  in  cylindrical  coordinates, 

2  r„,\2  r.Al        -■  ■     -2 


*  -* i  -d i  -r + m + im 


(-21) 


and 


(r^')2  =  t 


3J/ 
3r 


ail 


8z 


The  change  to  the  new  coordinates  n,  7.   is  quite  lengthy,  but  not 
difficult  otherwise.   The  results  finally  obtained  can  be  summarized  in 
the  following  form. 

r2    X  r    f    \ 

n  =  — r  G  (n) 
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(-22) 


(-23) 
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(-24) 


where  the  G's  are  known  functions  defined  as  follows. 
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Note  that  these  three  functions  are  all  finite  on  the  axis  n  =  fl 
and  that  they  all  vanish  at  infinity. 

To  find  the  averaging  distance  X      in  a  non-dimensional  form  appropriate 
to  this  problem,  we  evaluate  the  expression 


All2 


■  a,  ■  - 


/OO     A  .00      9 

a   rdr      /  G,(n)ndn 
0  '  0  1 


1   r(Ofi')2rdr   H   [C2(n)  +  C2(n)]r,dn 


(-28) 


2V '"        ^3 

It  should  be  observed  that  both  integrals  in  this  last  equation  are 
evaluated  from  r,  =  0  to  n  =  °°  over  the  surface  z  =  constant,  and  that 
in  this  case  the  quantity  z  cancels  identically  from  the  right  side  of 
the  equation.   There  is  no  integration  with  respect  to  z.   The  dimensionless 
averaging  distance  a.  is  constant  and  independent  of  z. 

The  integrals  involved  in  Eq  (-28)  can  easily  be  evaluated  numerically, 
thus  fixing  a...   It  is  found  that  a1  is  roughly  the  same  as  s,  the  jet 

width  parameter;  this  seems  reasonable. 

2 
A  similar  procedure  is  involved  in  determining  the  integrals  I   and 

2 
J  and  the  resulting  local  length  scale  parameter.   Also  the  effects  of 

polar  symmetry  must  be  taken  into  account  in  the  manner  of  Appendix  B. 

Again  the  integrations  are  from  n.  =  0  to  n  =  °°  over  the  surface  z  = 

constant,  and  again  the  variable  cancels  from  the  final  result.   Thus 
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Hence  the  dimensionless  loca]  length  scale  becomes 
\}2  2,  ,    1  fl2)    Il(n) 

z    =  a  (n)  =  ib  =  727; 

;  z  VJ  '   J.  (n) 

2         2 
The  integrals  I,(n)  and  J,(n)  must  of  course  he  evaluated  numerically 

It  turns  out  that  the  local  length  scale  a(n)  obtained  in  this  way  is  a 

slowly  varying  function.   The  value  remains  nearly  constant  at  ahout  0.1 

across  most  of  the  jet,  and  declines  only  slcvly  at  relatively  large 

distances  from  the  centerline. 

The  energy  equation  for  the  axi-symmetric  turbulent  jet  may  be 

written,  first  in  cylindrical  coordinates, 
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The  various  functions  vhich  appear  in  the  energy  eciuation  can  nov 
be  represented  as  functions,  of  the  similarity  coordinates  z  and  r,  as 
follows 
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Upon  substituting  the  expressions(-33)  into  the  energy  equation,  we 
obtain , 
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This  result  confirms  the  fact  that  all  terms  of  the  enerpy  equation 

-A 
are  proportional  to  z   and  that  the  above  energy  lav  is  indeed  compatible 

with  the  self-similarity  of  the  overall  solution.   In  order  to  expand 

this  so  as  to  be  able  to  eliminate  z  from  the  result,  we  male  use  of  the 

relation 
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In  this  way  we  obtain  the  general  result 
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Considerable  simplification  can  be  achieved  with  little  loss  of 
accuracy  by  neglecting  the  effect  of  the  slight  divergence  of  the  jet 
on  the  energy  equation.   This  amounts  to  the  assumption  that 


_3 
3z 
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(-37) 


In  this  case  the  energy  equation  simplifies  to 
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(-38) 


To  complete  the  solution  it  is  necessary  to  know  the  empirical 
functions  a,  6,  and  y   as  defined  by  Eqs.  (-18),  (-19),  and  (-20)  of 
Section  1.3.   In  this  case,  however,  there  arc  no  fixed  walls  involved, 
and  these  quantities  become  constants,  namely, 
a  =  0.065 


e=-3T7 


(-39) 
Y  =  1.4 
The  eddy  viscosity  t(n)  which  occurs  in  the  energy  equation  may  he 
found  in  either  of  two  ways.   Firstly,  we  may  use  the  solution  previously 
found  from  the  momentum  analysis,  by  numerical  integration  of  En.  (-15). 
Secondly,  we  may  use  the  basic  postulate  of  the  heuristic  model,  that  is, 


E  =  aA/2E  =  aa,/2L 

This  second  alternative  represents  a  somewhat  more  severe  test  of 
the  heuristic  model  and  was  adopted  for  the  numerical  example  here 
discussed.   Thus  substituting  Eq .  (-40)  into  (-38)  gives 
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The  energy  function  E  (n)  is  the  only  unknown  remaining  in  this 
equation.   The  boundary  conditions  are 

•3E. 


At    n  -  0 


3ri 


=  o 


(-42) 


At    n-*0  E+°°  (-43) 

Eq .  (-41)  must  be  solved  numerically.   It  happens  to  be  convenient 
for  this  purpose  to  use  the  substitution 

2EX  =  q*  (-44) 

The  equation  is  then  expanded  and  the  derivatives  replaced  by  appropriate 

finite  difference  approximations  in  a  straightforward  manner. 

The  boundary  condition  (-43)  at  infinity  must  be  satisfied  indirectly. 

Actually,  the  integration  starts  at  the  centerline  using  -r-1"  =  0  and 

^  '  o 

a  trial  value  of  q  .   A  value  of  q   is  found  by  iteration  which  yields 

o  o  J  J 

q»  0  at  some  large  value  of  n.   In  this  case  the  value  n  =  0.25  may  be 
regarded  as  "large.'   It  is  found  that  the  precise  position  of  the  outer 
boundary  makes  virtually  no  difference  in  the  results  as  long  as  it  is 
greater  than  about  0.20. 

The  energy  distribution  found  in  this  way  is  shown  by  the  solid  line 
in  Fig.  1.6-2.   The  agreement  with  the  experimental  results  of  Corrsin, 
Kef.  (7),  Laurence,  Ref.  (8),  is  seen  to  be  satisfactory,  especially  in 
view  of  the  degree  of  scatter  in  the  data  points. 

Once  the  energy  distribution  is  known  the  eddy  viscosity  is  fixed 
by  Eq.  (-40)  and  the  corresponding  shear  stresses  by  Eq .  (-20).   These  are 
shown  by  the  solid  line  in  Fig.  1.6-3.   The  dashed  line  shows  the  corres- 
ponding shear  stresses  as  computed  from  the  momentum  equation.   Data  points 
are  from  Corrsin,  Eef.  (9). 
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If  the  heuristic  model  were  completely  correct,  and  if  the  assumed 
Gaussian  velocity  profile  were  exact,  then  the  two  shear  stress  curves 
shown  in  Fig.  1.6-2  would  coincide  exactly.   As  it  is,  the  degree  of 
agreement  attained  is  considered  to  be  reasonably  satisfactory.   The 
discrepancy  between  the  two  curves  is  smaller  than  the  discrepancy  between 
the  experimental  points  and  either  of  the  curves. 

It  seems  probable  that  in  this  instance  the  theoretically  computed 
values  are  actually  more  accurate  than  the  experimentally  measured  ones; 
the  experimantal  measurement  is  innately  difficult  and  uncertain. 

It  may  be  concluded  that,  on  a  whole,  these  results  for  the  turbulent 
jet  substantiate  quite  well  the  proposed  heuristic  model  of  fluid  turbulence, 

1.7  CONCLUSIONS  AND  RECOMMENDATIONS 

It  is  concluded,  on  the  basis  of  the  evidence  available  so  far,  that 
the  proposed  heuristic  model  is  basically  adecmate  for  determining  the 
principle  flow  characteristics  in  the  general  case  of  inhomogeneous  and 
non-stationary  turbulence  in  incompressible  flow. 

It  is  recommended  that  the  present  model  be  applied  also  to  other 
cases  including  pipe  flow,  two  and  three  dimensional  wakes,  boundary  layers 
and  to  some  examples  of  unsteady  flow.   The  aim  would  be  to  refine  the 
present  unified  theory  and  extend  its  range  of  applicability. 

It  is  also  recommended  that  further  experimental  information  be 
obtained  in  connection  with  those  aspects  of  the  model  for  which  the 
present  data  are  insufficient.   These  aspects  include,  for  example, 
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(a)  the  relation  between  the  local  length  scale  X  of  the  mean  flow 
and  the  correlation  length  X*  of  the  turbulence, 

(b)  the  generalized  three  dimensional  stress/strain  rate  relations 
which  actually  exist  in  regions  of  strong  anisotropy  such  as 
in  the  flow  near  a  wall,  and 

(c)  the  influence  of  various  key  parameters  on  the  rate  of  dissipation 
of  turbulent  energy  into  heat. 
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Part  II 
Theory  and  Background  of  the  Heuristic  Model 

by 
T.  H.  Gawain,  D.Sc. 


A3 


2.1   THE  FUNDAMENTAL  PROBLEM 

Consider  a  general  turbulent  flow  field,  one  which  is  both  in  — 
homogeneous  and  non-stationary.   How  can  the  turbulent  structure  of 
such  a  field  be  adequately  analyzed  and  described? 

In  principle,  it  is  possible  to  formulate  a  valid  general  method 
for  computing  numerically  the  detailed  solution  u.   (x,t)  for  a  particular 
realization  of  the  field,  subject  to  specified  boundary  conditions  and 
initial  conditions.   Once  these  conditions  are  fully  specified,  the 
equations  of  continuity  and  motion  suffice  to  establish  a  fully  determinate 
and  detailed  solution.   In  theory,  this  solution  procedure  could  be  repeated 
very  many  times  to  generate  a  whole  ensemble  of  such  detailed  realizations. 
All  realizations  would  be  compatible  with  the  same  general  macroscopic 
boundary  conditions  and  initial  conditions.   However,  they  would  differ 
more  or  less  randomly  in  regard  to  the  fine  scale  details  of  the  boundary 
and  initial  conditions,  and  therefore  also  in  the  fine  scale  details  of  the 
resulting  motions. 

Once  such  an  ensemble  of  detailed  solutions  were  available,  any 
desired  statistical  features  of  the  turbulence  could  then  be  computed 
from  this  ensemble  in  a  completely  determinate  way. 

It  is  significant  that  the  foregoing  method  does,  in  principle, 
satisfy  the  mathematical  requirements  for  a  determinate  solution.   This 
solution  is  based  solely  and  strictly  on  the  general  equations  of  contin- 
uity and  motion  (the  Navier  Stokes  equations),  and  on  the  pertinent  bound- 
ary and  initial  conditions  for  the  particular  case  in  question. 
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In  practice,  however,  the  volume  of  detailed  calculations  required 
to  carry  out  the  above  hypothetical  scheme  of  calculation  for  any  case 
of  practical  interest  is  simply  overwhelming.   It  lies  vastly  beyond  the 
capabilities  of  any  computer  yet  known  or  contemplated.   The  basic  reason 
for  this,  of  course,  is  that  the  detailed  space/time  structure  of  turbulence 
is  extremely  fine  grained  and  enormously  intricate.   Hence  for  all  practical 
purposes,  the  above  hypothetical  solution  procedure  must  be  ruled  out. 

There  is  a  classical  method  for  attempting  to  circumvent  the  above 
difficulty.   Since  it  is  impractical  to  solve  the  detailed  equations  and 
average  the  results,  the  classical  method  averages  the  basic  equations 
themselves.   It  then  attempts  to  compute  the  desired  average  characteristics 
of  the  flov  from  these  averaged  equations.   Now,  if  the  basic  equations 
were  linear,  this  technique  would  indeed  solve  the  problem.   Unfortunately, 
they  are  not.   As  a  result  of  the  non-linearity,  the  averaging  process 
introduces  additional  unknowns,  the  Reynolds  stresses,  into  the  basic 
equations.   Consequently,  there  are  finally  more  unknowns  than  available 
equations.   Therefore,  the  averaged  equations  do  not  comprise  a  closed 
and  determinate  set.   This  constitutes  the  so-called  closure  problem  of 
turbulence  theory. 

It  should  be  emphasized  that  the  deficiency  which  arises  from 
averaging,  the  equations  is  fundamental.   No  amount  of  mere  manipulation 
of  the  equations,  however  involved,  can  remedy  this  defect.   The  averaging 
of  the  equations  leads  to  a  definite  loss  of  essential  information  in 
comparison  with  the  information  implicit  in  the  original  equations. 

Therefore,  so  long  as  we  are  constrained  to  work  only  with  the 
averaged  equations  instead  of  with  the  originals,  it  becomes  unavoidably 
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necessary  to  invoke  additional  hypotheses  in  order  to  compensate  for  the 
lost  information  and  to  reestablish  a  determinate  set  of  equations. 
Unfortunately,  there  is  no  known  V7ay  to  establish  the  required  hypotheses 
deductively  from  the  original  fundamental  equations  themselves.   Consequently, 
it  is  initially  far  from  clear  just  what  auxiliary  hypotheses  are  most 
nearly  adequate  and  satisfactory. 

Obviously,  the  basic  justification  for  any  particular  hypotheses 
finally  adopted  must  therefore  rest  on  a  comparison  of  the  consequences 
deduced  from  these  hypotheses  with  experimental  observations.   This  is 
the  accepted  scientific  method  of  inductive  inference.   While  it  is  to  be 
regretted  that  the  problem  apparently  cannot  be  solved  on  a  wholly 
deductive  basis,  there  are  no  a-priori  reasons  for  despairing  over  the 
prospects  for  ultimate  success  by  the  inductive  method. 

It  should  also  be  pointed  out  at  this  juncture  that  the  averaged 
equations  of  motion  describe  exactly  the  effect  of  the  Reynolds  stresses 
on  the  mean  motion.   However,  they  give  no  information  whatever  concerning 
the  reciprocal  effect  of  the  mean  motion  on  the  Reynolds  stresses.   Yet 
some  definite  reciprocal  effect  of  this  kind  unquestionably  does  exist. 
This,  therefore,  is  precisely  the  problem  with  which  our  auxiliary 
hypotheses  must  necessarily  deal. 

2.2  STATISTICAL  DESCRIPTION  OF  TURBULENCE 

Consider  an  arbitrary  and  general  flow  field  which  is  neither 
homogeneous  in  space  nor  stationary  in  time  in  regard  to  its  turbulent 
structure.   Let  u^   (x,t)  denote  any  one  particular  realization  of  this 
field.   This  means  that  the  velocity  components  u.   Cx, t)  represent  one 
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particular  detailed  solution  of  the  equations  of  continuity  and  motion; 
this  solution  is  compatible  with  the  specified  macroscopic  boundary 
conditions  and  initial  conditions.   There  exist,  of  course  an  infinity  of 
such  possible  individual  realizations  which  satisfy  these  general 
conditions.   The  separate  realizations  differ  from  one  another  only  with 
respect  to  their  fine  scale  details.   We  denote  each  realization  ty  a 
distinct  value  of  the  index  n,  which  is  taken  to  be  an  integer.   Thus  we 
may  have,  for  example,  n  ■  1,  2,  3,...N  where  N  is  some  indefinitely 
large  number.   Such  a  very  large  and  representative  collection  of  individual 
realizations  we  term  an  ensemble. 

Fortunately,  the  particular  details  of  any  one  individual  realization 
are  seldom  of  direct  interest.   What  is  usually  required  is  statistical 
information  concerning  the  principal  features  of  the  ensemble  as  a  whole. 
In  particular,  with  respect  to  any  variable  of  interest,  what  we  would 
generally  like  to  know  are  not  the  individual  values  for  particular 
realizations,  but  rather  the  following  items  of  information.   These  are 
listed  in  order  of  decreasing  urgency.   Thus 

(a)  The  average  value  of  the  variable  over  the  ensemble. 

(b)  The  statistical  variance  of  the  quantity  with  respect  to  the 
ensemble  average. 

(c)  The  detailed  probability  distribution  of  the  quantity  over 
the  ensemble. 

For  example,  if  U'    (x,t)  is  the  instantaneous  total  velocity  at 
space/time  point  x,t  for  realization  n,  we  require  the  corresponding 
ensemble  mean 


i(x,t)  -  Aim|  I     U^n;Cx,t)  (2.2-1 


N 

I 

N-*»  '  n=l 
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This  result  defines  the  mean  flow  field.  Of  course,  the  turbulent 
fluctuations  of  the  nth  realization  are  defined  as  the  deviations  from 
this  mean.   Thus 

u£n)(x,t)  =  l^(n)(x,t)  -  Ui(x,t)  (-2) 

From  this  definition,  it  follows  that  the  turbulent  fluctuation  u^   has 
zero  mean  over  the  ensemble. 

The  corresponding  variance  of  this  velocity  fluctuation  is  given  by 


~2 1    N    (n)2 

u.(x,t)  ■  Aim   , N_-.N  I     u.    (x,t)   (No  sum  on  index  i)      (-3) 

N-*»   ^    >    n=i 

Finally,  the  probability  function  which  defines  the  detailed  statistical 
dlsttributlon  of  this  random  variable  u.  over  the  ensemble  may  be  denoted 
by  the  functional  4'[u   (x,t)].   Here  the  quantity  fdu  represents  the 
probability  that  for  a  realization  n,  drawn  at  random  from  the  ensemble, 
the  velocity  component  u.  at  point  x,t  will  have  a  value  which  lies 
between  u  and  (u  +  du.).   If  the  space/time  coordinates  x,t  be  he 
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constant  at  prescribed  values,  f(u.)  reduces  to  a  corresponding  ordinary 

probability  density  distribution.   Hence  the  functional  ¥  amounts  to  an 

ensemble  of  probability  distributions,  one  for  each  space/time  point  of 

the  field. 

If  the  turbulence  happens  to  be  stationary  with  respect  to  time  in 

its  statistical  features,  then  it  becomes  permissible  and  convenient  to 

replace  the  ensemble  average  by  a  simple  time  average  over  a  single 

realization.   The  realization  index  n  may  therefore  be  dropped.   For 

example,  the  mean  velocity  distribution  may  now  be  redifined  as  follows. 

+T 
U  (X)  -  Aim  y=  /  U'(x,t)dt  (-4) 

T-*»    -T 
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Similar  changes  apply  to  the  definitions  of  variance  and  probability 
distribution. 

Likewise,  if  the  turbulence  is  homogeneous  with  respect  to  some 
line,  surface  or  volume,  the  ensemble  average  may  be  replaced  by  the 
corresponding  spatial  average. 

In  the  following  analysis,  the  average  will  always  be  an  ensemble 
average,  unless  specified  otherwise.   Also  the  discussion  will  center 
largely  on  certain  important  averages  and  variances,  but  will  not  deal 
specifically  with  the  corresponding  detailed  probability  distributions. 

2.3  CORRELATION  FUNCTIONS 

In  order  to  formulate  an  adequate  theory  of  turbulence,  it  is 
necessary  to  define  certain  basic  concepts  which  characterize  the 
turbulence  in  the  general  vicinity  of  any  prescribed  space/time  point. 
In  this  connection,  the  two  point  velocity  correlation  tensor  is  of 
basic  importance. 

Let  P  denote  a  reference  space/time  point  with  coordinates  x,t 
regarded  as  arbitrary  but  fixed.   Let  P1  and  P"  be  two  associated  points 
with  coordinates 

x'  =  x  +  Ax 
t1  =  t  +  At 
and  (2.3-1) 

x"  =  x  -  Ax 
t"  =  t  -  At 
The  quantities  2Ax,  2At  represent  the  space  and  time  separation  between 
the  points  P'  and  P"  and  are  regarded  as  variable.   Note,  however  that 

49 


P"  and  P"  are  always  disposed  symmetrically  equal  and  opposite  with 
respect  to  the  fixed  reference  point  P. 

Components  of  the  turbulent  velocity  fluctuations  at  the  variable 
points  P'  and  P"  are  denoted  by 

u!  =  u  (x  +  Ax,  t  +  At)   at  P' 
and  (-2) 

u"  =  u,(x  -  Ax,  t  -  At)   at  P" 
The  mean  product  of  these  two  velocity  components  is  the  tensor 


u'.u'.'  =  u.(x  +  Ax,  t  +  At)  u.  (x  -  Ax,  t  -  At)  (-3) 

where  the  overbar  denotes  an  ensemble  average. 

At  zero  separation,  the  above  tensor  reduces  (upon  reversal  of 
algebraic  sign)  to  the  ordinary  Reynolds  kinematic  stress  at  the  reference 
point  P.   Thus 


t..  =  -u.u.  =  -u,  (x,t)u  .  (x,  t)  (-4) 

!J      x  J     k      J 


By  reason  of  its  symmetry  the  stress  tensor  u.u,  =  u.u.  has  at  most 

i  j    J  i 


six  distinct  components.   The  general  velocity  correlation  tensor  u!u'.' , 

however,  is  not  necessarily  symmetrical  with  respect  to  an  interchange 

of  the  indices  and  may  therefore  have  up  to  nine  distinct  components. 

Upon  contracting  the  Reynolds  stress  tensor,  we  obtain  the 

physically  important  quantity 

2  _      ~~2   ~~ 2   ~" 2 

q  (x,t)  -  u.u.  =  u-  +  u?  +  u_  =  2E(j?,t)  (-5) 

here  E  represents  the  ensemble  mean  turbulent  energy  at  the  point  x,t  and 

2 
q   is  the  corresponding  mean  square  velocity. 

We  now  define  the  space/time  correlation  tensor  R. .  merely  by 


normalizing  u'u"  with  respect  to  the  value  of  q  at  the  fixed  reference 
point  P.  Thus 
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u!uV   u  (x  +  Ax,  t  +  At)u.(x  -  Ax,  t  -  At) 

P   =  P   (x,t;Ax,At)  =  -±rl   =  -i - J 

J     J  q  q  (x,t) 

(-6) 
It  is  of  interest  to  compare  the  above  definition  which  we  have 
adopted  with  certain  alternative  forms  which  are  sometimes  encountered. 
One  such  alternative  is 


u.(x,t)u.(x  +  Ax,  t  +  At) 
p.   =_i J (_7) 

1J    q(3T,t)q(x  +  Ax,  t  +  At) 


Another  alternative   is 


u.(x,t)u.(x  +  Ax,    t  +  At) 

R.  .    -    ,X  J  ,  =^        (No   sum  on   i  or  j)  (-8) 


y  u.(x,t)  /uT(x  +  Ax,    t  +  At) 


The  form  we  have  chosen  to  preferable  because  it  is  related  in  a 
simpler  fashion  to  the  energy  spectrum  of  turbulence,  and  because  it  has 
certain  useful  symmetry  properties  which  will  be  pointed  out  presently. 

Two  important  special  cases  of  the  foregoing  are  the  spatial  and 
temporal  correlation  tensors,  defined  respectively  as  follows: 


,  v                  u .  (x  +  Ax,  t)u.(x  -  Ax,  t) 
R{S)    -  R   (X,t;Ax,0)  --^ J 

2  J  q   C*,t)  (_9) 


m  u  (x,t  +  At)u  (x,t  -  At) 

R)\}    =  R..(x,t;0,Ax)  =  — = J 

ij     1J  q^(x,t) 

Of  particular  significance  are  the  first  invariants  of  these  two 
tensors,  defined  as  follows: 
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,  v  u  (X  +  AX,  t)u  (X  -  AX,  t) 

Rii  "  R±1Cx»t;fia,o)  -  — 2~^ 

q  <*•*>  (-10) 


/Ts  u.(x,t  +  At)u  (x,t  -  At) 

R)/    -  R   (x,t;0,At) 


ii     iiv  '  '  '   '  2,_  x 

q  (x,t) 

Of  course,  the  summation  convention  applies  as  usual  to  the  repeated 

(s) 
indices  in  the  foregoing  expression.  Note  that  R..   is  an  even  function  with 

(T) 
respect  to  the  independent  variable  Ax.  Likewise  R.  .   is  even  with  respect  to 

At.   These  important  symmetry  properties  are  a  consequence  of  defining 

R   according  to  Eq.  (-6),  rather  than  according  to  the  alternative  (-7) 

or  (-8). 

Because  of  their  fundamental  character  and  intimate  relation  to 
the  turbulent  energy,  the  scalar  correlation  functions  R;.  and  P.. 
provide  the  natural  bases  for  establishing  appropriate  length  and  time 
scales  of  the  turbulence  in  the  neighborhood  of  the  space/time  point  x,t. 

Specifically,  let  dv1  =  dxidxldx'  represent  an  infinitesimal 
element  of  volume,  and  dt1  an  infinitesimal  element  of  time  at  tbe  variable 

point  P' .  Also  let 

2        2        2 
Ax'Ax  =  Ax  Ax  =  (Ax-)   +  (Ax?)   +  (Ax_) 

represent  the  square  of  the  spatial  separation  PP ' . 

Then  the  desired  characteristic  length  A*  and  time  t*  may  be  defined 

as  follows : 

.              /|R(f  |(Ax-Ax)dv' 
X*Z  =  X*Z(x,t)  - *±7-z (-11) 

/|R<f |(At)2dt' 

T*^  =  T*Z(x,t)  =  ^JZX (-12) 

/i*2}idt' 
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The  integrals  in  (-11)  extend  over  all  space  at  a  fixed  time  t. 
The  integrals  in  (-12)  extend  over  all  time  at  a  fixed  spatial  point  x. 
Of  course,  these  integrals  are  finite,  despite  the  infinite  domains  of 
integration,  because  of  the  rapid  decay  in  the  correlation  functions  at 
large  separations. 

Note  that  the  integrands  contain  only  the  absolute  values  of  the 
correlation  functions;  this  ensures  that  every  non-zero  element  of  the 

correlation,  even  if  locally  negative  in  some  places,  makes  a  positive 

2       2 
contribution  to  the  integrals  which  define  A*   and  t*  . 

The  numerators  of  Eqs.  (-11)  and  (-12)  represent  the  second  moments 

of  the  respective  correlation  functions.   They  are  therefore  analogous 

2       2 
to  moments  of  inertia.   Also,  the  quantities  A*  and  t*  are  analogous  to 

variances  in  probability  theory.   It  is  well  known  that  such  second 

moments  attain  their  minimum  values  when  evaluated  with  respect  to  the 

centroidal  point  of  the  distribution  in  question.   In  the  present  instance 

the  correlation  functions  are  symmetrical,  with  their  centroids  at  the 

respective  origins  of  coordinates,  namely  at  Ax  =  0  and  at  At  =  0.   Hence 

2       2 
the  characteristic  quantities  A*   and  t*  represent  the  minimum  possible 

"variances'1  of  the  respective  correlation  functions.   These  properties 

make  A*  and  t*  particularly  appropriate  measures  of  the  length  and 

time  scales  of  turbulence  in  the  vicinity  of  the  arbitrary  space/time 

point  x,t. 


2.4  SPECTRUM  OF  TURBULENCE 

The  space/time  correlation  tensor  P..  may  be  transformed  into  the 
wave  number/frequency  domain  by  means  of  multi-dimensional  Fourier 
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integrals.   Let  k\    k2  K3  ^e  tne  components  of  the  spatial  wave  number  vec- 
tor k.   Let  03  be  temporal  wave  number,  that  is,  angular  frequency.   Also 

let  dv  =  dK-dK  Hk.  represent  an  infinitesimal  element  of  "volume"  in 
k     12   3 

wave  number  space.   Note  therefore  that 

<*Ax  ■  k.Ax.  =  k-Ax..  +  k_Ax_  +  K-Ax_  (2.4-1) 

We  may  now  define  the  Fourier  transform  pair 


4>  .(x,t;K,w)  = 


« Ju  J 


[  /-  ».a-  a*n  +i(<*Ax  4-  wAt)  ,  ,,.  ,      ,  _* 

R. . (x, t;Ax,At)e  dv  dt       (-2) 


R   (x,t;Ax,At)  =  — - 
1J  /2tt  J  J 


*  i—   *.-   \  ~i(<*Ax  +  wAt),   ,        ,  ON 
4> .  .  (x,t;<  ,o))e  dv  dw       (-3) 


The  integrations  in  Eq.  (-2)  extend  over  all  of  space/time  Ax, At 

with  parameters  x,t  held  constant.   The  integrations  in  (-3)  extend  over 

all  of  wave  space  k,w  with  parameters  x,t  again  held  constant. 

The  spectral  correlation  function  $ .  .  like  its  Fourier  transform 

ij 

R. .  is  a  nine  component  tensor.   Naturally  these  two  related  functions 

exhibit  corresponding  symmetry  properties. 

The  first  invarient  4>   like  its  Fourier  transform  R..  is  a  scalar 

ii  li 

and  has  a  corresponding  basic  significance.   It  defines  the  character  of 
the  energy  spectrum  in  the  vicinity  of  the  space  time  point  x,t. 
Incidentally,  it  can  be  shown  from  energy  considerations  that  * . ,  must  be 
non-negative  everywhere. 

It  may  be  seen  that  the  functions  R   and  <t>   are  exactly  equivalent 
in  the  information  they  provide  concerning  the  structure  of  the  turbulent 
field.   The  reason  is  simply  that  if  either  one  of  these  functions  be 
specified,  the  other  one  can  be  found  from  the  above  Fourier  integral 
relations. 
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It  is  clear  therefore  that  fundamental  scales  of  length  and  time, 
which  were  previously  developed  in  terms  of  R    and  R    can  equally 
well  be  defined  in  terms  of  the  corresponding  transforms  $j.   and  *i.  • 

Because  of  the  above  equivalence  between  the  correlation  functions 

R..  and  the  spectral  functions  * . .  it  is  immaterial  in  terms  of  which  of 
ij  ij 

these  our  turbulence  theory  is  developed.  We  have  chosen  to  emphasize 
the  correlation  aspects,  and  it  is  seen  that  there  is  no  loss  of  infor- 
mation involved  in  so  doing. 

It  should  be  noted  that  the  nine  spectral  functions  P .  .  (x, t;Ax,At) 
are  ostensibly  eight  dimensional,  the  independent  variables  being  x..  ,  x_, 
x_,  t,  Ax-  ,  Ax_,  Ax_,  At.   Similarly,  the  nine  spectral  functions  <J> 
(x,t;K,u))  are  also  eight  dimensional,  the  independent  variables  being 

X-i,  X_,  ^o»  ^  >  ^  1  '  *^9  '    \  *     ^  * 

In  addition  to  the  above  eight  independent  dimensions  which  are 
explicit,  there  is  an  additional  dimension  which  is  implicit.   Recall 
that  the  correlation  functions  are  themselves  averages  theoretically 
computable  from  a  large  ensemble  of  individual  flow  realizations.  We  may 
assign  to  each  realization  a  superscript  index  n  which  varies  from  one  to 
some  very  large  number.   Then  n  is  equivalent  to  an  additional  independent 
dimension.   For  if  we  wish  to  know  not  just  the  ensemble  averages  them- 
selves, but  also  something  about  their  statistical  distributions,  then  it 
becomes  necessary  to  take  into  account  the  possible  variations  of  the 
independent  realization  n  from  the  ensemble  mean. 

On  the  other  hand,  any  specialization  of  the  conditions  reduces 
the  number  of  independent  variables  which  must  be  considered.   Thus,  for 
stationary  turbulence,  R.  .  becomes  independent  of  t  (but  not  of  At)   For 
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completely  homogeneous  turbulence,  R   becomes  independent  of  x..  ,  x„ ,  x^. 
Also,  for  Isotropic  conditions,  spherical  symmetry  prevails,  and  R. 

ceases  to  depend  on  Ax.,  Ax9,  Ax_  individually,  but  becomes  a  function 

I  2        2        2 

only  of  their  resultant v (Ax.)   +  (Ax_)   +  (Ax_)  .   Furthermore,  in  this 

case  the  number  of  distinct  components  of  R  .  reduces  from  nine  to  two, 

depending  only  on  whether  the  indices  i,j  are  equal  or  unequal.   These 

are  vast  simplifications.   They  account  for  the  fact  that  homogeneous 

and  isotropic  turbulence  has  been  studied  so  much  more  extensively 

than  any  other  case. 

It  should  be  noted  that  the  fundamental  phenomena  of  turbulence  as 

represented  by  the  ensemble  of  functions  u   (x,t)have  five  independent 

variables,  namely,  x..  ,  x_,  x_,  t  and  n,  and  three  dependent  variables,  namely, 

u. ,  u_,  and  u».   It  is  a  curious  fact  that  when  we  seek  to  make  these 

phenomena  more  comprehensible  by  means  of  the  correlation  functions  R.., 

the  number  of  dependent  and  independent  variables,  instead  of  being 

sharply  reduced,  are  actually  increased!   Of  course,  the  functions  R .  . 

have  a  regular  and  symmetrical  structure,  they  vary  smoothly,  and  die 

away  rapidly  at  large  separations,  whereas  the  u.  vary  in  an  extremely 

complex  and  erratic  way.   It  is  this  qualitative  difference  that  makes  the 

correlation  functions  seem  so  much  more  comprehensible. 


2.5  RELATION  BETWEEN  TURBULENT  CORRELATION  SCALES  AND  MEAN  FLOW 

It  is  customary  and  in  many  ways  convenient  to  resolve  fluid 
turbulence  into  two  distinct  complementary  aspects  which  we  then  refer  to 
as  "the  mean  flow"  and  "the  turbulent  fluctuations",  respectively.   While 
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such  a  distinction  is  useful  for  certain  conceptual  and  computational 
purposes,  it  should  not  be  forgotten  that  these  two  aspects  are  inextri- 
cably coupled  within  a  single  unified  process.   In  this  connection  it  is 
pertinent  to  note  that  while  we  often  find  it  convenient  to  introduce 
this  distinction  between  mean  flow  and  turbulent  fluctuation  into  the  basic 
Navier  Stokes  equations,  the  distinction  is  not  innate  within  the  original 
equations  themselves.   This  is  shown  by  the  fact  that  these  equations, 
being  non-linear,  cannot  be  fully  decoupled  into  two  independent  sets, 
one  dealing  exclusively  with  the  mean  flow,  the  other  dealing  exclusively 
with  the  turbulent  fluctuations. 

The  foregoing  fact  provides  a  useful  hint  in  connection  with  the 
problem  of  formulating  appropriate  auxiliary  hypotheses,  as  required  for 
supplementing  the  averaged  equations  of  motion.  We  have  seen  that  any 
such  hypotheses  will  naturally  be  formulated  in  terms  of  the  fundamental 
physical  parameters  which  characterize  the  overall  turbulence  in  the 

vicinity  of  a  given  space/time  point  x, t.   Among  these  quantities  are 

2 
the  mean  square  velocity  q  (x,t)  and  the  macroscales  of  turbulence  in 

space  and  time,  as  expressed  by  the  previously  defined  parameters  X*  and 

T*. 

Of  course,  in  most  problems  of  practical  importance,  interest  centers 
chiefly  on  the  mean  flow  field.   One  reason  for  this  is  that  the  mean 
flow  is  usually  much  more  evident  physically,  much  more  convenient  to 
observe  and  measure,  seemingly  more  easy  to  understand.   We  tend  to 
regard  the  turbulent  fluctuations  as  a  nuisance  which  has  to  be  dealt  with 
in  order  to  get  to  the  mean  flow,  and  we  tend  to  confine  attention  to 
those  aspects  of  the  turbulence  whose  connection  to  the  mean  flow  is 
obvious  and  unavoidable. 
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More  often  than  not,  we  do  not  have  available  the  data  which  would 
be  needed  to  establish  such  largely  unseen  characteristics  of  the  turbul- 
ence as  the  spatial  macroscale  X*  and  temporal  macroscale  T*.   Therefore, 
it  is  of  doubtful  utility  to  formulate  any  heuristic  theory  of  turbulence 
directly  in  terms  of  these  quantities.   For  example,  what  good  does  it 
do  to  express  the  Reynolds  stresses,  say,  in  terms  of  X*  and  t*  if  the 
latter  quantities  are  initially  just  as  much  unknown  as  the  Reynold's 
stresses  themselves?   Consequently,  we  must  at  this  stage  of  our  knowledge 
replace  the  fundamental  macroscales  X*  and  t*  in  our  formulation  by  other 
quantities  which  are  more  or  less  equivalent,  but  which  are  far  more 
accessible  to  observation  and  calculation. 

Fortunately,  for  this  purpose  we  can  taVe  advantage  of  the  intimate 
connection  as  noted  above,  which  necessarily  exists  between  the  local 
turbulence  field  and  the  local  mean  flow  field.   Presumably  because  of 
this  connection,  these  local  scales  of  space  and  time  must  somehow  be 
reflected  in  the  local  characteristics  of  the  mean  flow,  and  must 
therefore  be  deducible  from  the  mean  flow. 

While  the  exact  nature  of  this  relationship  between  the  macroscales 
and  the  mean  flow  is  far  from  clear  we  can  legitimately  infer  that  the 
relationship  is  inherently  somewhat  local  in  character.   This  follows  from 
a  fact  universally  revealed  by  all  observed  correlation  curves.  These 
always  show  a  rapid  decrease  in  correlation  with  increasing  separation  in 
space  and  time.   Moreover,  the  correlations  become  negligible  at  separa- 
tions larger  than  about  two  or  three  multiples  of  X*  or  T*. 

Consequently,  we  infer,  or  at  least  we  hypothesize,  that  quantities 
more  or  less  equivalent  to  the  true  macro  seales  X*  and  t*  can  be  defined 
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and  computed  in  terms  of  the  observable  characteristics  of  the  local  mean 
flow.   Let  us  call  these  quantities  the  apparent  macroscales  and  designate 
them  by  the  symbols  X  and  t  to  distinguish  them  from  the  true  macroscales 
X*  and  t*  as  earlier  defined.   The  important  problem  of  just  hov  X  and  x 
should  be  defined  in  terms  of  the  observable  mean  flow  characteristics 
is  taken  up  in  a  subsequent  section.   For  the  present  we  observe  simply 
that  if  A  and  x  be  defined  in  some  suitable  manner,  then  it  seems  probable 
from  the  above  reasoning  that  the  dimensionless  ratios 

1*    -    £T  ('2> 

will  be  very  stable  and  well  behaved  universal  functions  which  vary  but 
slowly  if  at  all.   We  assume  that  for  many  practical  applications  it  will 
turn  out  to  be  premissible  to  approximate  these  universal  functions  by 
suitable  average  values,  and  to  treat  them  as  simple  quaii-constants.   The 
actual  numerical  values  of  these  ratios  is  regarded  as  a  matter  to  be 
determined  empirically.   Of  course,  these  numerical  values  will  depend  on 
the  exact  definitions  that  are  finally  adopted  for  X  and  x.   However,  the 
presumption  is  that  if  the  definitions  adopted  are  near  optimal,  the 
ratios  will  turn  out  to  be  roughly  of  order  unity. 

It  follows  from  the  foregoing  discussion  that  our  definitions  of  the 
apparent  macroscales  X(x,t),  x(x,t)  in  the  vicinity  of  an  arbitrary  space/ 
time  point  x,t  should  depend  only  on  mean  flow  conditions  in  a  finite 
region  surrounding  that  point.   This  principle  is  more  realistic  than  any 
method  which  seeks  to  define  these  macroscales  in  terms  of  mean  flow 
quantities  only  at  the  point  x,t  itself.   On  the  other  hand,  it  avoids 
the  opposite  extreme  which  would  hold  that  X(x,t)  and  x(x,t)  are  somehow 
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directly  dependent  on  conditions  at  all  points  in  the  flow  field,  however 
remote  in  space  or  time.   Our  proposed  model  accords  more  nearly  than 
either  of  these  alternatives  with  the  observed  fact  that  turbulent 
correlations  are  effectively  limited  to  certain  finite  separations. 

Naturally,  the  definitions  of   A(x,t)  and   i(x,t)  should  emphasize 
mean  flow  conditions  in  the  immediate  vicinity  of  the  point  x,t,  should 
give  progressively  less  weight  to  conditions  farther  away,  and  should 
finally  neglect  mean  flow  effects  at  points  sufficiently  remote  in  space 
or  time  from  point  x,t.   This  requirement  suggests  the  general  idea  of 
defining  A  and  t  in  terms  of  appropriate  weighted  averages  of  certain 
pertinent  mean  flow  quantities. 

The  preferred  choice  of  weighting  functions  for  this  purpose  would 
be  the  correlation  functions  R.,   and  R. .   themselves,  if  these  were 
known.   Of  course,  if  these  functions  were  known,  they  could  be  used  to 
find  the  true  macroscales  A*  and  t*  directly,  and  there  would  then  be  no 
need  for  weighting  functions.   However,  they  are  not  known.   Consequently, 
we  are  obliged  to  used  heuristic  substitutes,  preferably  functions  which 
resemble  the  correlation  functions  to  a  certain  extent.   Fortunately, 
great  accuracy  is  not  necessary  in  this  regard,  for  weighted  averages 
tend  to  be  relatively  insensitive  to  minor  variations  in  the  form  of  the 
weighting  functions  used. 

In  view  of  these  considerations,  we  choose  the  Gaussian  curve  as 
an  appropriate  form  of  heuristic  weighting  function.   It  expresses  the 
general  trends  of  correlation  in  a  suitable  manner,  and  has  various 
convenient  and  known  mathematical  properties  as  well.   Consequently,  we 
finally  write  the  space/time  weighting  function  as  follows 
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Ax'Ax  _  (At)2 
~     2       2 
W  =  e     A       T  (-3) 

The  reason  for  inserting  the  numerical  factor  3  in  the  first  term  of  the 

exponent  will  be  explained  presently. 

Of  course,  for  purely  spatial  weighting,  we  may  delete  the  temporal 

term  in  the  exponent  and  set 

_  Ax'Ax 

W(s)  =  e     *2  (-4) 

Similarly,  for  purely  temporal  weighting,  we  may  delete  the  spatial  term 

in  the  exponent  and  set 

(At)2 


W(T)  =  e    ^  (-5) 

The  above  weighting  functions  W  '   and  W    may  be  regarded  as 

(s) 
initial  heuristic  approximations  to  the  true  correlation  functions  P, . 

li 

(T) 
and  K    ,  respectively.   Note  that  all  four  of  these  functions  equal 

unity  at  zero  separation. 

If  we  now  replace  the  true  correlation  functions  by  their  heuristic 

Gaussian  approximations  in  Eqs.  (2.3-11)  and  (2.3-12),  we  obtain  the 

interesting  results 

2    2 
X*Z  =  \*  (-6) 

2    2 

T*/  =  T  (-7) 

This  shows  that  for  Gaussian  correlation  functions  with  exponential 
constants  A  and   x,  the  true  macroscales  X*  and  x*  may  be  simply  identified 
with  these  exponential  constants.   Of  course,  this  simple  identification 
is  possible  only  if  the  numerical  factor  3  is  inserted  into  the  spatial 
exponent  as  indicated.   This  arises  from  the  fact  that  the  spatial  integra- 
tion of  Eq .  (2.3-11)  is  three  dimensional,  whereas  the  temporal  integration 
of  Eq.  (2.3-12)  is  only  one  dimensional. 
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Recall,  however,  that  the  actual  correlation  functions  are  not 
necessarily  Gaussian  and  are  in  fact  unknowns  in  the  present  analysis. 
Moreover,  for  the  reasons  previously  given,  the  quantities  X  and  t  will 
be  defined  not  in  terms  of  the  correlation  functions,  but  heuristically 
in  terms  of  the  mean  flow.   In  fact  the  basic  reason  for  introducing  the 
above  weighting  functions  in  the  first  place  is  that  they  are  needed  in 
connection  with  these  heuristic  definitions  of  X  and  t.   However,  the 
weighting  functions  in  turn  contain  X  and  t  as  parameters.   It  is  there- 
fore apparent  that  the  actual  calculation  of  A  and  t  will  necessarily 
involve  some  type  of  iteration  process. 

For  certain  purposes  in  the  subsequent  analysis  it  will  prove 
convenient  to  utilize  a  normalized  form  of  the  above  Gaussian  weighting 
function.   In  this  connection,  consider  the  following  integrals,  namely, 

K  =  //Wdv'dt'  (-8) 

H(s)  =  //W(s)dv'dt'  (-9) 


and 


H(T)  =  yvw(T)dv'dt  (-10) 


As  usual,  the  integrands  extend  over  the  entire  space/time  domain.   Of 
course,  the  integrals  themselves  are  always  finite  even  if  the  domain  of 
integration  extends  to  infinity. 

We  now  define  the  normalized  weighting  functions  as  follows: 


w 
v  "h 

(s) 
w         ■ 

w(8> 

H(s) 

(T) 
w          = 

(T) 

1! 

(-11) 

(-12) 
(-13) 
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When  normalized  in  this  way,  the  above  weighting  functions  have  the 
following  convenient  properties,  that  is, 

//wdv'dt1  =  1  (-1A) 

//w(s)dv*dt'  =  1  (-15) 

//w^dv'dt'  =  1  (-16) 


2.6  LOCAL  LENGTH  AND  TIME  SCALES  OF  THE  MEAN  FLOW 

We  have  been  considering  any  arbitrary  and  general  turbulent  flow 
field  of  which  the  mean  flow  may  be  either  steady  or  unsteady.   With  every 
space/time  point  of  this  field  may  be  associated  two  parameters  X  and  t 
which  express  respectively  a  length  scale  and  a  time  scale  characteristic 
of  the  mean  flow  pattern  in  the  vicinity  of  the  point.   Our  problem  is  to 
devise  suitable  explicit  definitions  of  these  parameters  in  a  way  which 
meets  various  essential  requirements.   Some  of  these  have  already  been 
discussed.   There  are  some  additional  requirements  as  well,  among  which 
are  the  following: 

(1)  The  parameters  X  and  t  must  be  true  scalars,  and  therefore  invarient 
with  respect  to  any  rotation,  reflection,  translation  or  acceleration 
of  the  reference  axes. 

(2)  These  parameters  should  also  be  everywhere  continuous,  finite  and 
positive  (except  possibly  in  certain  limited  singular  regions,  such 

as  at  a  solid  boundary) . 

2      2 

(3)  If  feasible,  the  mean  flow  quantities  X   and  x   should  preferably  be 

2       2 
related  to  the  correlation  quantities  X*  and  t*  ,  respectively,  such 

2       2 

that  the  ratios  — r  and  ——   are  close  to  unitv,  or  at  least  so  that 

X*2     t*2 


these  ratios  remain  as  nearly  constant  as  possible, 
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These  requirements  can  be  satisfied  as  follows.   Let 
Ul  =  ui(x1»x2'x3't)  =  ^(^O 
represent  the  velocity  field  of  the  mean  flow.   Then  the  corresponding 
strain  rate  tensor  is 

'3U.    31" 


(2.6-1) 


ij 


3x 


3x . 

r 


(-2) 


Of  course,  we  are  considering  here  only  the  incompressible  case, 
for  which  the  first  invariant  vanishes,  that  is 


31' 

rii  =  2iax 


-   2 


rv3U.i 


3x. 


3U, 


9x. 


3U. 


3x 


3J-* 


0 


(-3) 


Consequently,  there  is  no  net  dilatation,  and  the  Y.,    represent  purely 
distortional  effects. 

We  now  define  a  generalized  strain  rate  ft,  and  a  generalized  strain 
rate  gradient  ft'  by  means  of  the  following  expressions,  namely 

fi2  ■  iVu  <-*> 


ft 


,2   (3ft 


3x. 

l 


3ft 


3x 


(-5) 


From  these  definitions,  the  following  relation  can  be  obtained.   It 
will  prove  useful  in  the  subsequent  work.   Thus 


«*'>2--iF 


A  3x 


3ft 


3x 


(-6) 


The  factor  1/2  which  appears  in  Eq .  (-4)  calls  for  a  brief  explanation. 
Consider  a  simple  plane  parallel  shear  flow  such  as  illustrated  in  Fig. 
2.6-1.   Let  the  velocity  be  parallel  to  axis  x  ,  and  let  the  shear  gradient 
be  in  the  direction  of  x  .   Then  r   =  T        =  V        =  Tl^   =  V        =  0.   Now 
the  only  non-zero  shear  strain  rate  is  T  .  =  r_..  .   Then  applying  Eq .  (-A) 


6A 


U(y) 


///////////////////////// 


In  this  case  the  definitions  of  il  and 
£V  reduce  to 


/// 


ii  =U=  (dU/dy) 
£i'=U"Md2U/dy2) 


FIGURE  2.6-1     EXAMPLE:  SHEAR  QUANTITIES  XI  AND  ii1 

IN   PLANE    PARALLEL   FLOW 
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we  find  that,  because  of  the  presence  of  the  factor  1/2,  we  obtain  the 

simple  result 

ft  ■  r 
"    12 

In  other  words  the  presence  of  the  optional  factor  1/2  identifies  ft 

with  the  ordinary  shearing  rate  In  a  simple  plane  shear  flow.   Eq.  (-4) 

may  be  regarded  as  an  extension  of  this  comcept  to  a  general  state  of 

strain  rate. 

It  can  be  shown  that  while  U.  and  F, .  depend  on  the  orientation  of 

the  reference  axes,  the  quantity  ft  is  invarient  in  this  regard,  and  is 

12    2 
therefore  a  true  scalar.   The  same  is  true  of  the  quantities  ft  ,  ft  ,  ft' 

2  2      2 

and  (ftft1)  .   Hence  X     and  i   can  be  conveniently  defined  in  terms  of  any 

independent  pair  of  these  quantities. 

To  illustrate  this  idea,  consider  the  expressions 

2    2  ft2 
A   =  C  H_  (_8) 

A  ft'" 
2    2  1 

2      2 
where  C,  and  V.*"   are  dimensionless  constants.   The  combinations  shown  may 

X  T 

be  seen  to  define  local  scales  of  length  and  time,  X   and  t.   In  fact, 
Eq.  (-8)  will  be  recognized  as  the  expression  for  the  von  Karman  mixing 
length,  or  rather,  as  a  generalization  of  the  concept,  for  von  Karman's 
original  definition  was  restricted  to  the  special  case  of  parallel  flow. 

One  difficulty  with  the  above  approach,  however,  is  that  the  functions 

2      2 
X  and  t  so  defined  tend  to  fluctuate  too  wildly  instead  of  being  smooth 

2      2 
and  well  behaved.   Moreover,  at  certain  singular  points  where  ft  or  ft1 

2     2 
are  locally  zero,  the  quantities  X     or  t  may  reach  values  of  zero  or 

infinity,  neither  of  which  is  premissible  or  physically  meaningful. 
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The  cure  for  this  difficulty,  of  course,  is  to  replace  the  local 

2       2 
point  values  of  Q  and  ft'  by  suitable  space/time  averages  evaluated  in 

the  vicinity  of  the  point.   Thus  we  may  write 

Q.2   =  J7wfi2dv'dt'  (-10) 

U12   =  J7wft  ,2dv'dt'  (-11) 

wbereupon 

y2   -   r2  9'2  (    -i-n 

A     -  C     — =^  (-12) 

X   ft'2 

2  1 

t     -  C     =4  (-13) 

n 

2      2 
Ihe  new  values  of  A   and  x   obtained  in  this  way  will  now  vary 

smoothly  as  required,  and  will  never  be  either  zero  or  infinite.   The 

2      2 
dimensionless  constants  C,  and  C  remain  to  be  specified. 

A      t  ' 

Note  that  the  normalized  form  of  the  weighting  function  is  used  in 
the  integrals  (-10)  and  (-11) .   This  is  important  because  it  ensures  that 
the  integrals  remain  relatively  insensitive  to  the  precise  form  of  the 
weighting  function. 

For  the  reasons  given  earlier,  a  Gaussian  curve  has  been  chosen  for 
the  weighting  function.   For  consistency,  the  characteristic  exponential 
constants  A  and  t  of  this  function  should  be  respectively  equal  to  the 
A  and  t  computed  from  Eqs.  (-12)  and  (-13).   In  other  words,  the 
averaging  intervals  in  length  and  time  should  be  approximately  equal  to 
the  characteristic  local  scales  of  length  and  time.   Some  iteration  will 
be  involved  in  satisfying  this  requirement. 

It  turns  out  to  be  useful  to  generalize  the  above  idea  even  further. 
For  this  purpose  we  define  the  following  integrals,  namely, 


67 


i2  =  ;/w(Q2)a(fi,2)bdv,dt'  (-14) 

J2  =  //vCfiVcfi'Vdv'dt*  (-15) 

where  the  exponents  a,  b,  c,  d  may  be  assigned  arbitrarily.   However, 
for  reasons  which  will  soon  be  evident,  we  stipulate  that 


D  = 


b       d 
(a  +  b)   (c  +  d) 


¥   0  (-16) 


Moreover,  to  avoid  improper  integrals  we  require  all  of  the  above  exponents 
to  be  non-negative. 

Now,  depending  only  on  how  the  arbitrary  exponents  a,  b,  c,  d  are 

2      2 
specified,  the  integrals  I   and  J  may  be  assigned  various  physical 

dimensions.   Consequently,  it  is  always  possible  to  define  a  mean  square 

2  2 

length  A   and  mean  square  time  t  by  expressions  of  the  following  general 

form,  namely, 

1  ?       ?    k  ?    $ 

\l  =  czxa  )  v  r  (-17) 

t2  =  C2(I2)m(J2)n  (-18) 

2      2 
where  C,  and  C  are  arbitrary  dimensionless  constants.   It  is  easily  shown 

A       T 

on  dimensional  grounds  that  the  four  exponents  k,  1,  m,  n  are  given  by 
the  following  formulas: 

k  =  -  £(c  +  d)  m  =  +  | 


(-19) 


£  =  +  ^(a  +  b)  n  =  -  - 

It  is  now  seen  that  Eq .  (-16)  represents  a  necessary  condition  to 
ensure  the  existence  of  the  above  solution. 
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The  foregoing  formulas  show  that,  subject  only  to  the  condition  (-16), 

2      2 
specification  of  the  six  arbitrary  constants  a,  b,  c,  d ,  C.  and  C  uniquely 

A  T 

2      2 
defines  the  quantities  A   and  t  .   These  six  ouantities  are  here  regarded 

as  empirical  constants  whose  values  can  in  principle  be  chosen  so  as  to 

2      2 
confer  upon  the  defined  quantities  X   and  t   a  wide  range  of  desired 

properties . 

Ideally,  it  would  be  desirable  to  choose  these  six  constants  in 

2      2  2 

such  a  way  as  to  make  X   and  t   agree  as  nearly  as  possible  with  X*     and 

2 
t*  ,  respectively,  over  a  broad  range  of  conditions.   In  principle,  this 

can  always  be  done  if  and  when  sufficient  experimental  data  becomes 

available  for  this  purpose.   At  the  present  time,  the  available  data  are 

insufficient  for  such  an  analysis.   It  therefore  becomes  necessary,  at 

least  for  the  time  being,  to  resort  to  expedients  of  a  more  heuristic 

nature,  as  will  be  explained  presently. 

2      2  2 

The  calculation  of  the  space/time  integrals  I   and  J  which  fix  X 

2 
and  t   is  seen  to  be  inherently  lengthy  and  complicated.   One  reason  is 

that  the  space/ time  integrals  involved  are,  in  their  present  form,  four 

dimensional.   The  second  reason  is  that  these  same  quantities  appear  not 

only  as  the  final  unknowns,  but  also  as  the  exponential  constants  in  the 

Gaussian  weighting  function.   Hence  our  next  endeavor  is  to  simplify  the 

analysis  as  much  as  possible  in  both  these  respects. 

In  this  connection,  we  observe  first  that  in  computing  the  integrals 

2      2 
I   and  J  ,  if  the  mean  flow  happens  to  be  steady,  the  integration  over 

time  cancels  identically  from  the  result,  and  the  integrals  reduce  to 

ordinary  volume  integrals.   Secondly,  we  notice  that  if  the  flow  field  also 
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happens  to  be  either  plane  or  axi-symmetric,  all  volume  integrals  can  be 
reduced  to  corresponding  surface  integrals.   Also  if  the  flow  field  happens 
to  be  of  the  type  known  as  self-similar,  all  volume  integrals  again 
reduce  to  surface  integrals.   If  the  flow  field  is  both  axi-symmetric  and 
self-similar,  then  a  double  reduction  is  possible,  and  the  volume  integrals 
can  be  reduced  to  line  integrals.   These  are  tremendous  simplifications. 

By  a  self -similar  flow  field,  we  mean  a  field  in  which  there  exists 
a  family  of  surfaces  across  which  the  mean  velocity  distribution  has  a 
common  pattern;  the  velocity  distribution  across  the  various  surfaces 
differs,  at  most,  only  in  the  respective  scales  of  length  and  time,  or 
what  amounts  to  the  same  thing,  of  length  and  velocity.   Many  common  flows 
in  jets,  wakes,  ducts  and  boundary  layers  are  known  to  be  self-similar. 
The  characteristic  surfaces  in  such  cases  are  conveniently  taken  as  cross- 
sections  approximately  normal  to  the  principal  flow  direction. 

This  consideration  of  self-similarity  shows  us,  thirdly,  that  for 
such  flows  each  characteristic  surface  possesses  its  own  overall  length 
and  time  scales;  these  are  constant  over,  any  one  such  surface,  but  may 
differ  from  one  surface  to  the  next.   In  addition  to  these  constant 
overall  scales  of  length  and  time  which  characterize  the  cross-section  as 
a  whole,  there  are  also  distinct  local  scales  of  length  and  time,  as 
considered  earlier,  in  the  vicinity  of  any  reference  point  on  the  surface. 
These  local  scales,  unlike  the  overall  scales,  vary  smoothly  from  point  to 
point  along  the  cross-sectional  surface. 

Of  course,  even  when  we  have  taken  these  factors  into  account  to 

2      2 
simplify  our  formulas,  the  definitions  of  A  and  t  will  remain  indeter- 
minate until  the  four  exponents  a,  b,  c,  d  and  the  two  proportionality 
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2      2 
constants  C  and  C  are  definitely  specified.   Obviously,  there  exists  an 

A  T 

enormous  amount  of  flexibility  here,  and  tbis  was  done  intentionally. 

The  flexibility  has  been  deliberately  provided  in  order  to  permit  matching 

2      2  2       2 

X      and  t   as  closely  as  possible  to  X*     and  t*  over  a  broad  range  of 

conditions. 

Naturally,  to  do  this  properly  will  require  much  data  beyond  what  is 
now  available.   Hence  at  the  present  time  our  scheme  has  far  more 
flexibility  than  we  really  need  or  can  properly  exploit. 

Consequently,  we  have  found  it  necessary  to  simplify  at  this  time  in 
what  may  seem  to  be  a  somewhat  arbitrary  manner.   Actually,  such  simplifi- 
cations were  guided  by  a  fair  amount  of  numerical  experimentation.   Accord- 
ing to  present  evidence,  they  seem  to  give  reasonably  satisfactory  agree- 
ment with  physical  experiment.   No  doubt,  there  remains  room  for  further 
improvement,  however. 

Specifically,  for  our  present  turbulence  model,  we  bave  cbosen  the 
following  numerical  values  for  the  four  arbitrary  exponents,  namely, 

a  =  2      b  =  0      c  =  1     d  -  1  (-20) 

Also,  it  proves  to  be  most  convenient  to  choose  the  proportionality 
factor  simply  as 

C*  =  1  (-21) 

2 
The  other  proportionality  factor  C   is  not  required  for  the  present  simpli- 
fied analysis. 

Let  us  summarize  the  situation  as  it  now  stands  after  the  foregoing 
simplifications.   The  integrals  below  extend  over  the  entire  flow  field. 
Thus 

Ax  =  (x'  -  x)  =  separation  variable  (-25) 
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W(x,x")  =  e       =  weighting  function  (-26) 

H(x)  ■  /Wdv1  =  normalizing  integral  (-27) 

w(x\x)  *  —  W  ■  normalized  weighting  function  (-28) 

n 

2  _>       A 
I  (x)  =  /wfi  dv1  ■  first  characteristic  integral  (-29) 

2  2 

J  (x)  =  /w(flfi  )  dv'  =  second  characteristic  integral  (-30) 

2.    I2 
X  (x)  =  — r-  =  local  length  scale  parameter  (-31) 

J 

Notice  that  Eq.  (-26)  ahove  differs  from  Eq.  (2.5-4)  in  that  the 
unessential  factor  3  has  been  eliminated  from  the  exponent  for  convenience. 
This  is  arbitrary  but  permissible;  it  amounts  to  a  rescaling  of  X. 

The  foregoing  relations  suffice  in  principle  to  define  the  required 
scale  factor  X  at  every  point  in  the  flow  field  as  required,  but  do  involve 
a  computational  difficulty.   The  problem  is  that  the  value  of  X  is  needed 
for  Eq .  (-26)  before  a  computed  value  is  available  from.  Ea .  (-31).   Hence, 
it  is  necessary  to  resort  to  an  iteration  process  such  as  the  following. 
Let  the  (n  +  1) th  approximation  to  X  be  defined  as: 

where 

£fl«  '  ^11  spaceV*'3'^'^'  <"33) 

Jn+l(i?)  "  ^11  spaceV^'>(^'>12<'v'  (-3*) 

The  weighting  function  W  is  defined  as  follows.   Let 

W  (*.*•)  =  exp  -  (*'  -  *)•(*'  -*>  (_35) 

Xn(x) 
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Then 


W  (x,x) 

w  (x,x)  =  — (-36) 

/  ,,       W  (x,x')dv' 
J  all  space  n 

We  now  define  (in  principle) 

A2(x)  =  Him  X2(x)  (-37) 

Advantage  should  be  taken,  of  course,  of  any  knowledge  concerning  the 

flow  to  hasten  the  convergence  of  the  X  's.   In  particular,  if  the  flow  is 

n 

self-similar,  such  as  the  flow  in  a  pipe,  duct,  or  boundary  layer,  or  the 

far  downstream  regior  of  a  turbulent  jet  or  wake ,  calculations  need  only 

be  performed  on  a  single  representative  cross-section  of  the  flow,  and 

the  results  rescaled  according  to  the  appropriate  nondimensionalizing 

parameters . 

As  an  initial  guess  to  start  the  iteration  process,  we  take  A   =  °° 

o 

everywhere.   Then  in  computing  A  ,  we  find  that  the  weighting  function 
W  is  equal  to  unity  everywhere,  and  so  we  obtain  simply 

A2  =    a11  SPace —  (-38) 

/  ..       (fifi')  civ 
•'all  space 

a  constant  independent  of  position.   Thus,  x9  is  the  first  nonconstant 
approximation  to  A  that  we  obtain. 

Numerical  experiments  with  various  flow  fields  have  shown  that,  in 
general,  the  convergence  of  the  A  's  is  extremely  rapid.   This  is  illustrated 
in  Table  2.6-1  which  shows  successive  numerical  approximations  to  the  A 
distribution  in  laminar  flow  in  a  two  dimensional  channel.   Consequently, 
for  the  purposes  of  the  calculations  presented  in  this  paper,  the 
approximation  was  made: 

A(x)  »  A2(x)  (-39) 
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Some  improvement  could  possibly  be  gained  by  continuing  to  higher 
approximations,  but  in  view  of  the  gross  simplifications  already  inherent 
in  the  model  as  a  whole,  the  authors  feel  that  such  refinement  is  not 
warranted.   Moreover,  absolute  accuracy  in  this  regard  is  not  essential. 
It  is  necessary,  however,  that  the  chosen  order  of  approximation  be 
adhered  to  consistently.   It  is  also  important  to  adjust  the  initially 
undetermined  numerical  constants  of  the  model  so  that  the  chosen  order  of 
approximation  gives  satisfactory  agreement  with  experimental  data.   The 
significant  numerical  constants  are  those  that  occur  in  the  heuristic 
expressions  for  the  quanities  a,  3,  and  y. 


TABLE  2.6-1 


Convergence  of  Macroscale  Distribution  in  a  Two-Dimensional  Channel 
with  Parabolic  (Laminar)  Velocity  Profile 


y/b  A  /b 

n  =  0  n  =  1  n  =  2  n=3  n  =  20 
0.0 
0.2 
0.4 
0.6 
0.8 
1.0 


2b  =  channel  width       y  =  distance  from  centerline. 
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0.775 

0.695 

0.674 

0.664 

0.775 

0.706 

0.693 

0.689 

0.775 

0.731 

0.727 

0.727 

0.775 

0.758 

0.758 

0.758 

0.775 

0.782 

0.781 

0.781 

0.775 

0.802 

0.800 

0.800 

2.7  REYNOLDS  STRESSES  AND  GENERALIZED  EDDY  VISCOSITY  CONCEPT 

In  dealing  with  an  incompressible  fluid,  it  is  customary  to  divide 
the  equations  of  motion  through  by  the  fluid  density.   As  a  result,  the 
viscosity,  pressure  and  stress  terms  all  appear  in  these  equations  divided 
by  the  density.   It  is  convenient  to  avoid  writing  out  the  density  explicitly 
each  time.   Therefore  we  define  the  quotient  of  viscosity  divided  by 
density  as  the  kinematic  viscosity,  and  assign  it  the  symbol  v.   Similarly, 
the  quotient  of  pressure  divided  by  density  becomes  kinematic  pressure, 
symbol  p.   Stress  divided  Vy  density  becomes  kinematic  stress,  symbol  t , . . 
Note  that  kinematic  pressure  and  kinematic  stress  have  dimensions  of 
velocity  squared. 

The  kinematic  Reynolds  stress  tensor  may  now  be  written 

t..  -  t.  -  -uTCT  (2.7-1) 

The  contraction  obtained  by  equating  indices  and  summing  over  the 
repeated  index  defines  a  scalar  invarient,  namely, 

Tij  =  Tll  +  T22  +  T33  =  _(U1  +  U2  +  U3}  =  "q2  =  ~2E  ("2) 

Now  the  Reynolds  stress  tensor  may  be  rewritten  as  the  sum  of  a 
purely  dilatational  tensor  plus  a  purely  distortional  tensor,  that  is, 

t..  -  "  iq26..  +  xl.  (-3) 

Conversely,  the  distortional  part  is  defined  by 

t!.   =  x..  +iq26..  (-4) 

Of  course,  the  individual  components  of  x!.  depend  on  the  orientation 
of  the  coordinate  exes.   However,  the  double  summation 

■|t'  t!.   =  t'2  (-5) 

2  ij  ij 
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is  invarient  with  respect  to  any  change  in  the  orientation  of  the  reference 
axes  at  any  point,  and  is  therefore  a  true  scalar.  We  may  think  of  t' 
as  representing  a  kind  of  generalized  scalar  amplitude  of  the  tensor  t'  . 

Of  course,  the  optional  factor  1/2  was  inserted  into  Eq .  (-5)  to 
maintain  consistency  with  the  corresponding  formula  for  the  strain  rate 
invarient  ft.   The  physical  justification  is  as  follows.   Consider  a  simple 
state  of  pure  shear  stress  such  that  x'   =  x'   =  xl»  =  x'   =  0.   Let  the 
simple  shear  x'   =  x'.  be  the  only  non-zero  stress.   In  this  case  Eq . 
(-5)  gives 

t*  =  t|2  (-6) 

Hence  the  invarient  x'  may  be  identified  with  the  magnitude  of  the  simple 
shear  stress  itself.  Eq .  (-5)  represents  an  extension  of  this  concept  to 
the  general  state  of  stress. 

A  somewhat  similar  analysis  may  be  made  of  the  rate  of  strain  tensor 
of  the  mean  flow  field.   Recall  the  definition 


•3U.    SIM 
l— J.  +  —1 


ii    I  3x ,    3x , 

For  an  incompressible  fluid,  the  associated  dilational  tensor  vanishes, 
that  is, 


(-7) 


Fii  =  Fll  +  r22  +  r33  =  2 


3U,    3U„    3U. 

=  0  (-8) 


fdV.         3U     9U. 


^3x-    3x_    3x«j 


Consequently  the  tensor  T..    is  itself  wholly  distortional.   An  amplitude 
scalar  ft  may  be  defined  through  the  invariant  relation 

\  Vij '  "2  ("9> 
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Next  we  wish  to  formulate  an  appropriate  way  to  express  the  general 

relation  that  may  exist  between  the  distortional  stress  tensor  x'   and  the 

distortional  strain  rate  tensor  r . . .   In  this  connection,  it  is  instructive 

ij 

to  consider  as  a  preliminary  model,  the  relation  which  is  known  to  exist 
between  the  ensemble  averaged  distortional  viscous  stress  (t'  )  and  the 
distortional  strain  rates  r . .  of  the  mean  flow.   This  is  given  by 


where  v  is  the  ordinary  molecular  kinematic  viscosity. 

This  model  suggests  that  we  attempt  to  express  the  distortional 
Reynolds  stresses  in  the  form 

t! .  =  eT. .  (-11) 

where  v  is  now  replaced  by  e,  the  so-called  eddy  kinematic  viscosity. 
However,  this  form  of  equation  implies  that 


x'    t'     t1     t*    t'    t' 
11    22    33    23    31    12 


(-12) 


r    r    r    r    r    r 

11    22    33    23    31    12 
Unfortunately,  there  is  no  a-priori  reason  to  believe  that  the  tensors 
t!.  and  r,.  will  necessarily  satisfy  the  five  similarity  constraints 
expressed  by  Eqs.  (-12).   We  therefore  require  a  more  general  relationship, 
one  which  remains  applicable  even  if  t!.  and  Y..    are  both  specified 
arbitratily. 

One  method  which  is  sometimes  suggested  for  dealing  with  this 
requirement  is  to  postulate  a  relation  of  the  form 

t'   =  e-.-.r,.  (-13) 

ij     ljkl  kl 

Now  the  eddy  viscosity  becomes  a  fourth  order  tensor.   Since  t'   and  T.. 

each  have  six  distince  components,  the  tensor  e    ...    has  thirty-six 


77 


components!   By  employing  tensor  arguments,  the  number  of  independent 

eddy  viscosity  components  can  actually  be  reduced  to  ten.   Furthermore,  if 

the  tensors  tI.  and  r..  are  arbitrarilv  specified,  Eqs.  (-13)  can  in  fact 
ij      ij 

be  reduced  to  an  equivalent  set  of  five  independent  equations  in  ten 
unknowns.   However,  these  equations  do  not  constitute  a  determinate  set. 
Clearly,  this  type  of  formulation  is  both  inadequate  and  overly  complex 
as  veil . 

A   somewhat  simpler  tensor  relation  has  been  suggested  by  liirt, 
Ref.  (A).   In  our  notation  it  becomes 

Ti:  *  i[£ikrkj  +  Vki1  ('U) 

While  this  is  an  improvement  over  Eq.  (-13),  it  is  not  free  of  difficulties. 

If  e,   is  taken  as  an  unsymmetrical  tensor,  Fq .  (-14)  fails  to  produce  a 

determinate  solution  for  e.  n  when  x!.  and  V..    are  specified.   On  the  other 

kl      ij      ij 

hand,  if  e...  is  taken  as  symmetrical,  the  relation  (-14)  loses  generality; 

there  are  then  certain  theoretical  combinations  of  t'   and  T.,  which  are 

ij      ij 

not  reconcilable  with  (-14)  no  matter  how  the  e.,  are  chosen. 

The  foregoing  difficulties  can  be  overcome  by  adopting  the  following 

relation  between  t'   and  T..,    namely, 

ij  ij 

t'     =  e[r,.  +  fif ..]  (-15) 

ij  ij  ij 

It  will   be  shown    that    the  quantities   c   and   f..    suffice   to  define   the 

ij 

relation  between  distortional  stress  and  strain  rate  in  a  fully  general 

way.   If  tI.  and  T        are  specified  arbitrarily,  a  corresponding  unique 

solution  can  be  found  for  c   and  f  ,  . .   Conversely,  if  e,  f..,  and  T..   are 

ij  ij       ij 

specified,  then  t'..  is  uniquely  determined. 
Solving  Eq.  (-15)  for  f..  gives 

f.,  =^[-t!.  -  r..]  (-16) 

ij   Q  Le  ij    ij 
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It  is  immediately  apparent  from  this  that 

tl}   =  ti±  (-17) 

and 

£ii  =  fll  +  f22  +  f33  "  °  <"18> 

However,  as  it  stands,  Eq .  (-16)  does  not  uniquely  fix  f..  hecause 
it  contains  the  as  yet  undetermined  parameter  e.   It  will  now  be  shown 
that  for  any  arbitrarily  prescribed  t!   and  r ,.,  a  corresponding  "preferred" 
value  of  e  can  be  uniquely  defined.   For  this  purpose  let 

t" .  =  efif  .  .  =  x ' .  -  er.  .  (-19) 

Here  t". ,  represents  the  discrepancy  between  t!.  and  eT . . .   We  wish  to 
ij  ij       ij 

choose  e  so  as  to  minimize  this  discrepancy  in  some  appropriate  sense. 
Hence  we  form  the  invarient  of  the  discrepancy  namely, 

t'.'.t'.  .  =  2t"2  =  [t!,  -  eI\.][t'   -  eV..] 
ij  ij  ij     ij   ij     ij 

=  t*  .t!  .  -  2eT  .  ,t!  .  +  e2r.  ,r. . 
ij  ij      iJ  iJ      ij  iJ 


or 


(-20) 

t"2  =  t'2  -  er..x!.  +  eV 
ij  ij 

2 
Now  to  minimize  the  invarient  t   of  the  discrepancy,  we  set 

2 

=  0  -  T..T,  .  +  2eQ2  =  0  (-21) 


3e         ij  ij 
whereupon 

F  =  13  iJ  =  ij  ij  =  ij  1J  r-22) 

0o2   r..r,.   r..r..  v  } 

2^     ij  ij    ij  ij 

This  simple  result  is  important.   It  is  the  basic  definition  of 
the  scalar  eddy  viscosity  e  for  the  case  of  general  stress  and  strain  rate. 

Now  multiplying  (-19)  be  V    .    and  using  (-22)  we  obtain  the  significant 

consequence  that 

r..f4.  -  0  (-23) 

ij  iJ 
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This  is  the  mathematical  expression  of  the  concept  that  the  tensors  T 

ij 


and 


f   are  'orthogonal".   Hence  our  basic  stress/strain  rate  law,  Eq 


(-15) ,  expresses  the  idea  that  the  Reynolds  stress  tensor  t ' .  consists  of 

a  tensor  eT    .   which  is  simply  proportional  to  T...   plus  a  tensor  eftf  . 
ij  '  ij  v  ij 

which  is  orthogonal  to  T  .  . . 

The  tensor  f..  has  nine  possible  components.   The  symmetry  constraint 

(-17)  reduces  this  to  six  distinct  components.   The  non-divergence  constraint 

(-18)  reduces  it  further  to  five  independent  components.   Tinally,  the 

orthogonality  constraint  (-23)  reduces  the  degrees  of  freedom  of  f . .  to  just 

four.   Hence  e  and  f..  together  possess  five  degrees  of  freedom.   This  is 

precisely  the  necessary  and  sufficient  number  to  express  the  general 

relation  between  distortional  stress  and  mean  flow  strain  rate.   If  we  add 

2 
7*—  as  a  sixth  degree  of  freedom,  this  then  allows  us  to  define  the  general 

relation  btween  overall  Reynolds  stress  and  mean  flow  strain  rate. 

Since  the  foregoing  equations  are  pure  tensor  relations,  they  remain 
valid  for  all  possible  orientations  of  the  coordinate  axes.   It  is  particul- 
arly convenient  to  choose  the  principal  axes  of  strain  rate  in  order  to 
bring  out  the  essential  concepts  in  the  simplest  possible  form.   We  use 
asterisks  to  designate  these  axes.   Also,  there  is  no  loss  in  generality 
in  ordering  the  subscripts  so  that  T*     £  r*  2.  T*.     Of  course  V*     =  T*^   = 
T*  =  0.   Moreover,  it  happens  to  be  useful  to  introduce  the  following, 
auxiliary  terminology.   Let 

rli-  A a  sln(f +  e* 

r*2  =  t»  B  sin(0  +  e*)  (-24) 

r$3  -  jl  0  sin(-  |S  +  8* 
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The  merit  of  Eqs.  (-24)  lies  in  the  fact  that  they  express  the  three 
quantities  r*  ,  r*  ,  F*   in  terms  of  just  two  parameters,  namely,  0.   and  G* 
Furthermore,  as  the  reader  may  easily  verify,  eqs.  (-24)  satisfy  identically 
the  following  important  constraints, 


r*  +  r*  +  r*  =  o 

11    22    33 


ir*2+  r*2   r,2 1    2 
2  !11    22    33   U 


(-25) 


The  parameter  6*  is  then  fixed  from  the  relation 


7  3  r*2 


ull   ]33; 


=  tan  6* 


(-26) 


It  is  equally  useful  to  adopt  a  similar  scheme  for  the  f* . .   Thus 


2   .   .  [2tt  ,  ,. 
fll  =v3  g  Sin  3~  +  * 


f22  =  7*3  g*  Sin(°  +  ^ 


(-27) 


i-j.     2      .  I      2tt    ,  . 
f33  =73g*  Bin-T+** 


whereupon 


11    22   r33 


(-28) 


12     2 

—  f*  +  f*  +  f* 

2   11    22    33 


=  g*' 


and 


(£fi  -  f§3> 


=  tan  ii"' 


(-29) 
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Now  the  orthogonality  constraint,  en.  (-23),  reduces  in  this  case  to 
the  form 

r*   f*  +  r*  f*  +  r*  f*  =  o  e-io") 

11  11    22  22    33  33  V    ; 

Upon  substituting  Eos.  (-24)  and  (-27)  into  (-30)  and  reducing,  we 
obtain  the  beautifully  simple  result 

^*  =  e*  ±  ~  (-31) 

This  too  expresses  a  kind  of  orthogonality  between  the  r.   and  the  f... 
With  ty*   now  known,  the  components  f*  ,  f$n>  f§^  are  all  fixed  as 
soon  as  g*  is  specified.   The  off-diagonal  terms  f£o>  f  t-i  >  anc*  ^io  may  ^e 
specified  individually.   The  final  stress/strain  rate  equations  for  these 
particular  axes  now  taVen  the  form 

tii  =  -ui  =  --h2  +  e[rii+  fiffi] 

T22  =  "U2  =  "  3  °2  +  G[r22  +  ff22] 

~2  12 

t*  =  _u  =  -  —  a     +   rr*  +  of*  1 

T33     3     3  q     l  33      33J 
T23  =  "V^  =  °  +  e[°  +  Cf23] 
T31  =  "V3!  =  °  +  ef°  +  fif*31l 

T12  =  "VI2  =  °  +  e[°  +  "f31] 

Note  that  the  principal  axes  of  stress  and  strain  will  coincide  if, 
and  only  if,  f*  =  f*x  =  f*2  =  0. 

The  further  reduction  of  these  equations  for  the  important  special 
case  of  a  plane  mean  flow  is  interesting.   Now  the  axes  x.  and  x_  are  in  the 
plane  of  the  flow,  and  x«  is  normal  to  this  plane.   Hence  we  have 


)*  =  0 


rii =  +r< 


(-33) 
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r*  =  o 

22 


r*  =  -r 
j33   i. 


(-33) 


Consequently  from  Bq .  (-31)  we  find  that 

V*  =  ±  |  (-34) 

There  is  some  choice  concerning  the  manner  of  handling  algebraic 
signs.   We  choose  to  drop  the  negative  sign  in  (-3^)  and  instead,  to  a] low 
g*  in  the  following  expressions  to  tale  on  positive  or  negative  values  as 
required.   Then 

f*  =  — ~   o*  f*   =  n 

rll     n   g  23 

f22  =  + j4  ?*  f31  =  arbltrary  (_35) 

f*   =  -  -^f  e*  f*   =  0 

r33     /I  g  12 

Of  course,  f*   =  f*   =  0  by  reason  of  symmetry.   Kence  g*  and  f* 
remain  as  the  two  arbitrary  degrees  of  freedom  which  define  the  coefficients 
f*  . 

The  final  stress/strain  rate  equations  become 

~2     12    nM     1   .. 
T*l  =  ~U1  =  "  I  q      [   3  R  ^ 

— 2     1   2  ? 

T22  =  ~U2  =  "  3  °     +   eC[°  +~1   g*] 

*     "2     12     .      1  (-36) 

T33  =  "U3  =  -  3  °     +  Efi[-1  ""3  g*] 


T*  =  -u  u   =  0  +  eTf* 
31     3  1         31 


T*   =  T*   =  o 
23    31 


Hence  the  overall  stress/strain  rate  relation  is  now  fully  specified 

1   2 
by  the  four  quantities  -r  q  ,  e ,  g*,  and  f*  . 
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In  the  case  of  axi-symmetric  mean  flows,  the  results  are  similar  to 
those  far  plane  flov,  except  that  6*  is  not  in  general  zero;  therefore  ty* 
must  be  found  from  Eq .  (-31)  and  the  f*   from  Eq .  (-27).   However,  we  still 

have  f*   =  f*  =  0  for  this  case. 

]   2 
The  basic  stress/strain  parameters  are  now--  q  ,  e,  g*,  f *  ,  and  ij.*. 

Note  from  Eqs.  (-20  and  (-31)  that  i|/*  is  a  known  function  of  the  strain 

tensor,  but  that  the  other  four  parameters  may  be  arbitrary. 

For  purposes  of  approximation  in  connection  with  our  heuristic 
turbulence  model,  we  have  decided  to  make  the  simplifying  assumption  f..  =  0. 
Actually  this  is  not  as  drastic  as  it  seems.   In  the  first  place,  as  we 
have  seen,  f..  embodies  at  most  only  four  degrees  of  freedom,  not  six  or 
nine.   Secondly,  we  are  mainly  concerned  with  plane  or  axisymmetric  mean 
flows,  for  this  f..  has  only  two  degrees  of  freedom,  namely,  g*  and  f*  . 
Thus  this  assumption  amounts  to  neglecting  but  two  degrees  of  freedom. 

Secondly,  an  analysis  of  the  momentum  transfer  involved  in  generating 

the  Reynolds  stresses  suggests  that  the  term  cT . .    accounts  for  most  of  the 

effect,  and  that  the  deviation  efif . .  should  indeed  be  small.  While  the 

ij 

experimental  data  which  is  available  is  not  decisive  on  this  point,  it  does 
tend  to  support  the  above  hypothesis.   It  does  seem  that  this  simplification 
may  become  somewhat  inadequate  in  the  immediate  vicinity  of  a  fixed  wall, 
but  it  appears  to  be  reasonably  satisfactory  elsewhere  in  the  flow. 

Under  these  assumptions,  the  principal  axes  of  stress  and  strain 
coincide,  and  the  equations  for  these  particular  axes  reduce  to  the  simple 

form 

*      1  2 
t*  -  -  -z  a     +   efi 

T*2  -  -  \   q2  +  0 

(-37) 
1*3  -  -  3  q   -  efi 

T*    =  T*   =  T*   -  0 
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The  foregoing  stress/strain  rate  law,  despite  its  apparent  simplicity, 
can  be  useful  for  calculating  the  Reynolds  stresses  which  correspond  to 
any  specified  mean  flow  strain  rates  only  if  the  value  of  e  can  somehow  be 
estimated  with  acceptable  accuracy.   Whether  this  can  be  done  is  an 
empirical  question,  capable  of  being  answered  by  suitable  experimental 
information.   The  experimental  data  so  far  available  is  not  yet  adeauate 
to  answer  this  need  fully.   Further  experiment  and  analysis  is  needed. 
Tentatively,  however,  we  are  justified  in  postulating  that  e  is  some  function 
of  the  characteristic  local  turbulence  parameters  q,  X,  and  t  as  previsouly 
defined,  perhaps  also  of  molecular  kinematic  viscosity  v,  distance  y  from 
the  nearest  wall,  and  possibly  of  some  as  yet  unknown  additional  parameters 
as  well.   This  can  be  written  symbolically  in  the  form 

e  =  e(q,  X,    t,  v,  y, )  (-38) 

Actually,  we  would  be  on  slightly  more  fundamental  ground  to  choose 
X*   and  t*  in  place  of  X   and  t.   Unfortunately,  X*   and  t*  are  not  usually 
known,  whereas  A  and  t  have  been  related  in  a  specified  manner  to  the  mean 
flow.   Consequently ,  for  pragmatic  reasons  we  replace  X*   and  t*  by  X   and  t 
in  the  above  expression. 

The  foregoing  empirical  law  which  governs  the  eddy  viscosity  e  can  be 

greatly  simplified  by  dimensional  analysis.   We  recal]  that  in  mechanics 

generally,  the  number  of  fundamental  dimesnions  is  three.   These  are  often 

taken  as  mass,  length  and  time.   However,  in  the  case  of  incompressible 

flow,  dividion  of  the  equations  of  motion  by  density  has  the  effect  of 

eliminating  mass  from  the  problem.   Our  final  equations  are  therefore  all 

kinematical  in  character,  with  length  and  time  as  the  only  fundamental 

dimensions.   It  follows  therefore  that  any  two  of  the  above  parameters  may 

be  chosen  as  the  necessary  and  sufficient  basis  for  the  non-dimensionalizing 

process . 
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From  a  purely  logical  standpoint,  it  would  seem  natural  and  convenient 
to  use  the  basic  local  scales  of  length  and  time,  A  and  T,  for  this  purpose. 
However,  there  are  grounds  for  believing  that  a   and  A  play  roles  which  are 
physically  more  important  than  that  played  by  t.   This  seems  to  be  especially 
true  for  mean  flows  which  are  steady,  or  which  vary  only  slowly.   Conse- 
quently, we  choose  q  and  A  and  our  non-dimensionalizing  parameters. 

As  a  result  of  the  foregoing  considerations,  the  empirical  eddy 
viscosity  law  may  now  be  rewritten  in  the  dimensionless  form 


-  a  =  a 


in  sA  y_ 

\    '   v  '   X'' 


(-39) 


Thus,  our  problem  has  been  simplified  to  that  of  finding  the  dimen- 
sionless coefficient  a  which  presumably  varies  as  a  function  of  the  dimen- 


sionless parameters  \-r~\    and 


qA 


and 


1 


,  and  possibly  also  of  one  or  two 


other  dimensionless  variables  as  well. 

In  practice,  the  usefulness  of  the  above  formulation  hinges  on  whether 
the  dimensionless  coefficient  a  does  in  fact  behave  in  a  reasonably  simple 
and  stable  manner.   Our  present  indications  are  that  this  is  indeed  so. 
Owing  to  limitations,  of  time  and  data,  we  have  not  fully  established  the 


CT 


and 


J 


ciA 


We  have 


actual  dependence,  if  any,  of  a  on  the  quantities 
found  heuristically,  however,  that  satisfactory  overall  results  can 
apparently  be  found  by  assuming  that  a  depends  only  on  the  wall  distance 
parameter  (y/A) .   The  following  expression  seems  to  give  reasonably  good 
curve  fits,  namely, 

A 


{1+e-l>   TJ 


(-40) 

Of  course,  this  expression  is  tentative.  More  data  is  needed  on  this 
point  in  any  case. 
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Naturally,  this  particular  expression  for  a  is  applicable  only  so 
long  as  we  specifically  retain  our  present  definition  of  A.   Any  modifica- 
tion of  the  definition  of  X    clearly  calls  for  a  corresponding  adjustment 
in  a.   Note  hovever  that  the  dimensionless  quantity 


X*q 


aX 
A* 


(-A1) 


should  remain  substantially  invarient  with  respect  to  such  details. 

However,  even  if  a  be  accurately  known,  and  X    as  well,  the  actual 

eddy  viscosity  e  cannot  be  evaluated  until  q  also  is  established.   Bearinp 

1   2 
in  mind  that  the  turbulent  energy  is  E  =.  -z   o  ,  it  is  apparent  that  the 

determination  of  q  must  involve  an  energy  analysis  of  the  turbulence. 

Physically,  this  is  fundamental.   It  ties  togetber  the  turbulence,  the 

eddy  viscosity  and  the  mean  flow  into  a  consistent  overall  picture.   It 

also  provides  the  additional  equations  which  are  necessary  to  form  a 

determinate  set. 

2.8  HEURISTIC  ENERGY  ANALYSIS  OF  THE  TURBULENT  FIELD 


Starting  with  the  Navier  Stokes  equations  of  motion,  it  is  possible 

to  develop  a  corresponding  exact  energy  equation,  as  shown  in  detail  in 

Appendix  A.      The  following  result  is  obtained,  namely, 

I  II  III  IV 


(^ 


3t 


=  h,.r..    -  ^Y..Y..    -|— [U.E]    -|- 
2    ij    ij         2    ij    ij         8x        l  9x 


iL 


P  + 


u  .u  .^ 
J    -1 


V 


9x 


iL 


3E 

9x. 


CT 


il 


3x 


(2.8-1) 
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This  equation  governs  the  evolution  over  time  and  distribution  over 
space  of  the  ensemble  mean  turbulent  kinetic  energy,  that  is,  of 

E  =  -  q  =2  ujuj  =  2  (U1  +  U2  +  U3}  (_2) 

The  five  terms  on  the  right  side  of  Eq .  (-1)  have  the  following 
physical  significance: 
Term  I  -  The  rate  of  energy  input  into  turbulence  as  work  done  by  the  mean 

flow  against  the  resistance  of  the  Peynolds  stresses. 
Term  II  -  The  rate  of  irreversible  energy  dissipation  into  heat.   Note  that 


the  quantity  y..y..  is  positive  definite.   This  ouantity  accounts 
for  the  local  entropy  production. 
Term  III-The  rate  of  convective  transport  of  turbulent  energy  E  by  the  mean 

flow  U. . 

l 

Term  IV  -The  rate  of  net  transport  by  turbulent  diffusion  of  the  turbulent 

pressure  plus  kinetic  energy  fluctuations. 
Term  V  -  The  rate  of  net  transport  of  energy  by  molecular  diffusion.   Note 

that  the  Reynolds  stresses,  as  well  as  the  energy,  appear  in  this 

term. 
While  Eq .  (-1)  is  theoretically  exact,  it  contains  several  terms 
which  cannot  be  evaluated  exactly  without  detailed  knowledge  of  the 
turbulent  fluctuations.   Such  detailed  information  is  not  practically 
attainable.   It  is  therefore  necessary  to  approximate  these  terms  by  suit- 
able semi-empirical  expressions.   These  expressions  will  necessarily  contain 
empirical  coefficients  whose  values  can  in  principle  be  chosen  so  as  to 
optimize  agreement  with  available  experimental  data. 

We  have  shown  in  the  preceding  section  that  the  Reynolds  stresses 
can  always  be  expressed  with  complete  exactness  and  generality  in  the  form 
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-vTaT.   =  t  .  .  =  -  ~  q26.  .   +   e[T  .  .   +   nf  .  .]  (-3) 

i  J    ij     3  '  ij      ij     ij 

where 
and 

r..f ..  -  o  (-5) 

Consequently ,  term  I,  the  work  input  expression  in  the  energy  equation 
reduces  exactly  to  the  result 

1  2 

4r..r.  .    E   efi  (-6) 

2  ij  ij 

however,  as  we  have  seen  earlier,  e  must  here  te  approximated  hy  an 
expression  of  the  form 

e  =  aAq  (-7) 

where  a  is  presumably  a  well  behaved  empirical  function,  the  expression 
for  which  has  been  given  earlier. 

Next  consider  the  dissipation  of  energy  into  heat,  as  given  by  term  II 

of  the  basic  energy  equation.   Dimensioral  considerations  suggest  that 

this  term  can  be  expressed  to  good  advantage  in  the  form 

,  2 

F'H  =  2^ijYij  =  -72  (-8) 

This  expression  amounts  to  a  definition  of  the  so-called  dissipation  length 
A  which  is  a  length  characteristic  of  the  microscale  of  turbulence.   This 
form  of  expression  would  be  very  convenient  if  the  dissipation  length  A 
happened  to  be  related  to  the  macroscopic  scale  A  in  some  fairly  simple  and 
invarient  way.   Unfortunately,  this  does  not  seem  to  he  the  case. 
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However,  experimental  studies  such  of  those  by  haufer,  Pef.  (*>)  , 

clearly  suggest  that  at  sufficiently  high  Reynolds  numbers,  the  heat 

dissipation  effects  tend  to  become  independent  of  Reynolds  number!   This 

implies  that,  contrary  to  Eq .  (-8),  the  optimal  form  of  expression  for  E 

H 

is  one  which  does  not  contain  viscosity  v  explicitly.   This  wil]  be  the 
case,  in  fact  if  we  assume,  for  example,  that 


y2  -  ■•  (£ 


(-9) 


where  B'  is  now  hopefully  some  well  behaved  empirical  function.   For 
combining  these  last  two  results  gives 

3 


h-nT 


(-10) 


which  is  not  explicitly  dependent  on  v.   This  form  is  similar  to  that 
suggested  by  Prandtl  in  Pef.  (2). 

The  fact  that  the  viscous  dissipation  occurs  in  a  form  vhich  does  not 
explicitly  contain  the  viscosity  seems  rather  paradoxical  at  first  sight, 
and  calls  for  some  comment.   It  can  be  given  a  reasonable  interpretation  by 
considering  turbulence  energy  processes  from  the  spectral  point  of  view. 
It  is  well  known  that  the  work  input,  term  I,  occurs  mainly  at  the  long 
wave  length  end  of  the  spectrum,  and  that  the  dissipation  to  heat,  term  II, 
occurs  mainly  at  the  short  wave  length  end.   At  high  Reynolds  numbers, 
there  is  a  broad  range  of  wavelengths  in  between,  the  so-called  inertial 
range,  for  which  both  work  input  and  heat  dissipation  are  negligible.   Each 
wave  length  within  this  inertial  range  simply  receives  energy  from  the 
longer  wave  length  components  and  transmits  this  energy  to  the  shorter 
wave  length  components.   It  would  appear  that  the  rate  of  transmission  of 
energy  thru  the  spectrum  is  largely  controlled  by  this  process  in  the 
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inertial  range,  rather  than  by  the  viscous  dissipation  process  itself. 
Apparently  the  viscous  dissipation  rate  easily  adjusts  itself  as  may  be 
required  to  dissipate  all  of  the  energy  coming  through  the  inertial  range. 
Since  the  viscous  process  is  not  rate  controlling,  the  effect  of  viscosity 
tends  to  drop  out  of  the  experimental  picture. 

These  experimental  facts  and  theoretical  considerations  suggest 
that  the  heat  dissipation  should  be  analyzed  by  the  same  general  method 
previously  employed  in  connection  with  the  eddy  viscosity.   For  example, 
we  might  assume  tentatively  that 

EH  =  ER(q,  A,  t,  v,  y, )  (-11) 

Then  non-dimensionalizing  on  the  basis  of  q  and  X  as  common  parameters,  we 
obtain 

•AE, 


i . 


H 


3 

>  q  - 


-"■-•'fr^  *■••••)  <-12> 
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which  is  seen  to  be  equivalent  to  Eq.  (-10) 

This  result  suggests  that  t '  is  a  well  behaved  universal  function 
which  can  be  established  empirically.   However  sound  in  principle,  this 
idea  is  rather  difficult  to  apply  at  the  present  time  owing  to  the  paucity 
of  suitable  experimental  data.   It  does  seem  plausible  that  for  mean  flows 
which  are  not  too  strongly  unsteady,  the  effect  on  8'  of  the  parameter 

might  prove  to  be  negligible.   On  the  other  hand,  efforts  to  establish 

alone  have  not  as  yet  proved  conclusive  one  way 
or  the  other.   This  is  due  in  part  to  limitations  in  the  data  and  in  part 
to  limitations  in  time. 

Instead  of  following  the  above  line  of  investigation  to  the  end,  we 
discovered  the  following  more  promising  modification.   Numerical  experimenta- 
tion suggested  that  it  might  provide  a  more  satisfactory  alternative. 
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l>  ) 

3'  as  a  function  of  — ^- 

v 


Recall  that  the  quantities  A  and  t  in  our  model  are  themselves  computed 
from  the  integrals  I  and  J.   hence  it  amounts  merely  to  a  rearrangement  of 
form  to  rewrite  Eq.  (-11)  as  follows: 


eh  =  EH(q»  x»  J»  v>  y» ) 


(-13) 

Now,  non-dimensional izing  on  the  basis  of  q  and  J  as  common  parameters,  we 
obtain 


(q  J) 


=  e 


L(qJ) 


L_  Ifil!l1/3   fi." 
.2/3'  v^J  ]    '  I  2 


1/3 


(-14) 


Actually,  the  significance  of  the  dimensionless  group  on  the  left  was 

established  by  certain  somewhat  speculative  physical  arguments  which  will 

not  be  repeated  here.  (See  ref.  11)   It  was  found,  however,  on  the  basis  of 

numerical  investigation,  that  this  particular  dissipation  parameter  8  seems 

to  remain  more  nearly  constant,  and  to  agree  better  with  physical  data 

than  does  the  parameter  6*  resulting  from  the  earlier  analysis.  Hence 

this  is  the  form  that  was  retained  for  our  present  heuristic  turbulence 

model. 

Here  8  is  a  slowly  varying  dimensionless  function.   The  dependence 

5>,l/3 


of  B  on  the  quantities 


(qJ) 


-r~rx   and  — 
2/3     v 


-j—  ]    was  not  actually  investigated, 


due  to  time  limitations.   However,  it  was  found  heuristically  that  6  seems 


to  depend  only  on  the  relative  distance  v*r 


to  the  nearest  wall,  if  any, 


and  that  far  from  any  walls  0  becomes  essentially  constant. 

The  following  empirical  equation  seems  to  produce  reasonably  good 
curve  fits,  namely, 


?-i 


-1 


B  - 


3.7 


1  +  e 


(-15) 
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Some  data  which  tends  to  confirm  this  formulation  is  given  in 
Appendix  C.   In  any  case,  it  is  clear  that  the  energy  dissipation  to  heat 
needs  further  investigation. 

Next  consider  term  IV  of  the  basic  energy  equation.   It  is  useful 
for  purposes  of  discussion 


u  .u  ,-\ 
P  +-4" 


=  u.p  +  u  -J-l 


(-16) 


The  second  of  these  terms  clearly  represents  the  net  transport  of  the 


fluctuating,  turbulent  kinetic  energy 


ru  .u  .^ 


by  the  action  of  the  turbulent 


velocity  fluctuations  u  .   Clearly  the  turbulent  kinetic  energy  — ^ 


u  ,u  .^ 


is 


both  a  scalar  and  an  extensive  property.  It  may  be  regarded  as  something 
which  is  physically  transported  along  with  the  element  of  fluid  mass.  In 
this  respect  it  is  quite  analogous  to  other  extensive  properties,  such  as 
salinity,  for  example.  It  is  customary  and  appropriate  to  represent  such 
net  turbulent  transport  by  an  appropriate  empirical  diffusion  coefficient. 
Hence  we  write 


u.u.1       ->   fu.u  > 

3x .  I  2 
l v 


i  9E 


3x 


(-17) 


i' 


On  the  other  hand  the  pressure  velocity  correlation  term  u.p  is  in 
a  somewhat  different  category  owing  to  the  fact  that  pressure  p  is  an 
intensive  property.   It  is  not  an  extensive  property  which  can  be  regarded 
as  being  transported  along  with  the  fluid  mass.  We  note  however  that  p  is 
a  dependent  variable.   This  means  that  whenever  the  spatial  distribution 
of  the  velocity  fluctuations  u.  is  specified,  the  corresponding  spatial 
distribution  of  p  is  likewise  fixed.   Hence,  looking  at  the  matter  statisti- 
cally, we  can  say  that  p  is  correlated  in  some  fashion  with  u..   In  fact, 
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taking  into  account  dimensional  considerations,  we  can  improve  this  state- 


J_J. 
2 


For  example, 


ment  and  say  that  p  is  statistically  correlated  with 
in  the  limiting  case  of  an  inviscid  fluid,  Bernoulli's  equation  informs 
us  that  regions  of  higher  than  average  kinetic  energy  will  tend  to  be 
regions  of  lower  than  average  pressure,  and  vice  versa.   Consequently  we 
are  justified  in  writing  the  pressure  term  in  the  analogous  form 


.,  3E 


u.p  =  -e 

i        3x 


(-18) 


xJ 


Now  the  overall  effect  becomes 


»t(p  + 


U  .U  "\ 

J  1 


=  -U,r  +  e') 


3E 


3x 


(-19) 


normally  c '  is  a  positive  quantity.   However,  because  of  the  generally 
negative  correlation  between  pressure  and  kinetic  energy,  we  expect  e"  to 
be  usually  negative  and  smaller  than  e'.   Hence  (e'  +  e')  should  be 
positive  but  smaller  than  e'  alone.   Moreover,  e"  and  e'  are  both  of  the 
same  dimension  as  the  kinematic  eddy  viscosity  e.   Mso  e'  can  be  expected 
to  be  roughly  of  the  same  order  of  magnitude  as  e.   For  these  reasons,  it 
is  advantageous  to  write. 

(e,:  +  c')  =  ye  (-20) 

where  y   is  not  a  dimensionless  coefficient  roughly  of  order  of  magnitude 
unity. 

Hence  finally  we  write 

'3E 


u  ,U."\ 

J  J 


UilP+  2   J 


=  -ye 


ox 


(-21) 


Various  numerical  experiments  suggest  that  the  dimensionless  empirical 
coefficient  y»  like  the  coefficients  a  and  f,  depends  only  on  the  distance 
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y  to  the  nearest  fixed  boundary,  if  any.   The  following  expression  has  been 


found  to  give  reasonably  good  results 

>2 


f"1 


y   =  1.4  -  0.4e 


(-22) 


Turning  next  to  term  V  of  the  basic  energy  eauation,  the  molecular 
diffusion  effect,  it  may  be  said  that  this  will  normally  be  very  much 
smaller  than  the  other  terms  in  the  equation.   It  may  safely  be.  neglected 
everywhere  except  in  the  very  thin  laminar  sublayer  region  along  any  solid 
boundary.  We  shall  ignore  it  in  our  mode]. 

Finally  we  combine  the  foregoing  steps  to  obtain  an  overall  heuristic 
energy  equation  as  follows.  Let 


q  -  2E 
Then  the  energy  equation  becomes 

,v7/6Tl/3 


(-23) 


3t 


=  efi  -B(2E)  '  J" 


i  iL 


3E 


3x 


(-24) 


where  e,  6,  and  a  are  fixed  by  the  following  empirical  expressions,  namely, 

2 


A72~E 


=  a  =  0.065  <  1  =   e 


{■ 


f-1 


-'f-i 


2      -1 


3.7 


1  +  e 


} 


Y   =   1.4   -  0.4e 


-If" 


(-25) 


(-26) 


(-27) 
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APPENDIX  A 


Derivation  of  Exact  Energy  Equation  for  Turbulence 


Equation  of  notion 

Let  prime  marks  on  a  symbol  be  used  to  denote  the  sum  of  the  mean 
value  of  that  quantity  plus  the  turbulent  perturbation.   That  is,  ]et 


u!  -  u.  +  u. 

J   j   J 


Also  let 


P*  =  P  +  p 


T  1  .  =  T  .  .  +  t .  . 
ij    ij    ij 


denote  viscous  shear  stresses. 

The  equations  of  continuity  and  motion  now  become 


su: 

i 

3x . 


=  0 


(*-l) 


(-2) 


(-3) 


3T!  ,-v 


3x 


ox.   11 
l    J 


(-A) 


The  corresponding  energy  eouation  is  found  by  multiplying  the  equation 


of  motion  by  U'.   Thus 
J 


9 

rU'.US 

J  J 

3t 

I    2      j 

=  u! 


m 


u!  |—  (D!u!) 


j^8xi  j    j  3x."  i  j 

Now  we  substitute  Lqs.  (-1)  and  (-2)  into  (-3),  (-4)  and  (-5),  then 
expand,  average  and  simplify.   All  terms  linear  in  perturbation  quantities 
vanish  in  the  averaging.   Therefore,  the  averaged  equations  of  continuity, 
motion  and  energy  become,  respectively, 


(-5) 
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3U. 
1 

3x. 

1 


-  0 


(-6) 


3b 


3t 


31 


JJ. 


3 


3x.      3x.   it     it 
l  J  x    J      J 


(-7) 


.  (U.U.    u.u.^ 

3   J  1  +  J  J 
3t   2      2 


-  U. 


-3T 
U  i — il 


Jl3x, 


+  u 


3t.  .^ 


i  3x. 


l  y 


+  C.  | — (u.u.  +  uTTT) 

]  3x.   l  j    i  j 


i  3x.   2 
l 


U  .U  *\ 

31 


-  u.u  . 


3U 


i  i  3x. 
j  v   1 


fU.U.U  -\ 

1  3  i 


3x. 

l 


(-8) 


Also,  multiplying  the  averaged  equation  of  motion  (-7)  through  by  I 
gives  the  energy  equation  for  the  mean  flov,  namely, 


3t  2 


3T.  . 
U  I — il 


i  1  3x . 
j  v  x 


-  u.  ^ — (u.u.  +  uTu:) 

i  3x.   l  j     l  j 


(-9) 


Now  subtracting.  (-10)  from  (-9)  gives  the  energy  equation  for  the 
turbulence  itself,  that  is 


3t 


3t 

ii 


3x. 


r  u.u.n     r3uo 

3x. I  l   2   J    l  i  3x. 
xu  J  l   i; 


ru  .u  .u . 
'  x  J  J 


3x.l    2 
i 


(-10) 


It  is  convenient  to  symbolize  the  mean  turbulent  kinetic  energy,  the 
instantaneous  kinetic  energy,  and  the  Reynolds  stresses,  respectively,  by 


E  = 


u  .u  . 

=   J  J 

2 


u  .u 


E'  = 


(-11) 
(-12) 


T  .  .  =  -U.U . 
1J        1  J 


(-13) 


Therefore  the  energy  equation  becomes 


3E 

at 


=  u . 


3t 


_il 


1  3x . 

l 


~ (U.E)  +  x.  . 

3x.   l      ij 


3U 


3x. 

l 


x 


(-14) 
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The  viscous  term  can  now  be  developed  further  as  follows.   Consider 
the  identity 


5t 


11 


Jl3x. 


=1x7(uitij)  ~  'ijlBx.j 

1     J      M  r 


3u 


(-15) 


Now  the  viscous  stress  perturbation  may  be  written 


t .  .    =    -p 6  .  .    +    VY  .  . 

1J  1J         1J 


(-16) 


where 


Y«  = 


9Uj         3u, 

+ 


9x.         9x . 

1    i        r 


(-17) 


and 


3u. 

l 

3x. 


=   0 


(-IP) 


Therefore 


u . t . .   =  -u .  p  +  v 


3E 


3t.  . 
1J 


3x .         3x 


3    ' 


(-19) 


t.  . 


dU  .'i 


ii \ 3x. 


=   0  +   v   y. 


:-u 


iJt^iJ 


(-20) 


Substituting    (-19)    and    (-20)    into    (-15)    gives 


3t.  . 


j    3x . 

l 


3 

3x 


.  oE 

-u.p  +    v  - — 
ir  3x 


t  .  .\-\  /-3u  .^| 


3 


(-21) 


Substituting  this  back  into  (-14)  and  rearranging,  we  find 


3E 

_^t 


ii I 3x. 


3U.-J 


ij  I  3x. 


3u  .^ 


3x.   l 


(-22) 


_3 
3x~ 


—  u 
: .   : 


u  .u 


P  + 


:D 


3x 


l  _ 


3E 

3x. 

l 


m 
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Two  further  simplifications  are  advantageous.   Recall  that 

•31'.    31' 


iJ 


1 


3x^    3x. 
i     J 


Noting  further  that  x   and  T..  are  both  symmetrical  tensors,  we  write 


or 


r3U.        31! 
i 


T.  .T.  .    =    T,  t|-r-=  +  -r-= 


=    2l.  . 
1J 


^3U.a 
_1 


3x . 


rdV 


iJ 


3x 


=  ^x.J 


iJ         2   iJ    « 


By  an   exactly  analogous   argument  we   find    that 


(-23) 


(-24) 


(-25) 


'3u  .^ 


Y, 


■  Fij^ij 


Ij^x^ 
Consequently,  the  turbulence  energy  equation  becomes  finally, 


(-26) 


3E 
3t 


2  ij  ij 


T^lj 


+ 


3x 


1  L. 


3E 


^7<uiE> 


3t.  ,>-r 


3 
3x 


u.  |p  +  J»J 


iL 


3x 


3x. 
J 


(-27) 


It  is  important  to  notice  that  this  equation  is  complete  and  exact. 
Each  term  has  a  definite  physical  significance  as  discussed  in  the  main 
body  of  the  report.   The  method  of  evaluating  certain  of  the  terms  heuris- 
tically  is  also  discussed  in  the  main  report. 
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APPENDIX  B 
Local  Length  Scale  for  Axi-Symmetric  Near  Flows 

Consider  the  following  sequence  of  relations  pertaining  to  the 
heuristic  model,  where  the  various  integrals  are  summed  over  the  entire  flow 

field. 

Ax*  Ax 


9 

xl 

W  =  e        =  Weighting  function  (B-l) 

H  =  /Wdv'  =  Normalizing  integral  (-2) 

W 
w  =  77  =  Normalized  weighting  function  (-3) 

n 
9        / 

I  =  /wfl  dv'  =  First  characteristic  integral  (-4) 

2  2 

J  =  /w(f-fi')  dv'  -  Second  characteristic  integral  (-5) 

2   I2 
A   =  —  =  Local  scale  parameter  (-6) 

J 
For  axi-symmetric  flows,  the  volume  integrals  (-2),  (-4) ,  and  (-5) 
reduce  to  surface  integrals  over  the  meridional  plane. 

In  this  case  the  volume  element  may  be  written  in  the  form 
dv'  =  r'de'dr'dz' 
Also  in  cartesian  coordinates 

Ax-Ax  =  (x*  -  x)2  +  (y'  -  y)2  +  (z1  -  z) 2  (-7) 

We  may  express  the  coordinates  of  the  fixed  point  x  and  of  the 
variable  point  x '  in  terms  of  cylindrical  coordinates  as  follows. 
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x  =  0 
y  =  r 
z  =  z 


y'  =  r'  cos 


(-8) 


Upon  substituting  (-8)  into  (-7)  and  rearranging,  we  can  write 

Ax«Ax  =  £2  +  2rr*(l  -  cos  6*)  (-9) 

where 

C2  =  (r'  -  r)2  +  (z'  -  z)2  (-10) 

2 
Now  let  us  consider  the  typical  integral  I   as  defined  by  Eq .  (-4) . 

By  making  use  of  the  foregoing  relations,  this  integral  can  he  reduced  to 

the  required  form.   In  addition,  it  is  convenient  to  utilize  the  following 

auxiliary  variable,  namely, 


2  _  rr 


(-11) 


Eq .  (-4)  is  now  expressible  in  the  form 


(r    ^ 


2tt 


■  M  /ff  l  lJ  fe  J  .*< 


-  cos  e')dei 


/rWdr'dz'     (-12) 


Now  since  Q   is  a  function  of  z'  and  r1  only,  but  independent  of  6 '  , 
the  above  integration  with  respect  to  6'  can  be  made  once  and  for  all. 
This  leads  naturally  to  the  definition  of  the  auxiliary  function 

2tt 


♦(c)  «t4  I  e"2^  (1  "  cos  e,)de' 


(-13) 


The  function  tjj  defined  above  must  be  computed  numerically.   The 
factor  —  was  inserted  on  the  right  so  as  to  achieve  the  simplification 
that  at  large  values  of  z,   the  function  i    approaches  the  limit  unity,  that 
is,  that 
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Aim  i>   =  1 

It  then  happens  that  for  small  values  of  Z, 
£im  i|j  =  2/¥  z, 


(-14) 


(-15) 


The  function  i>   is  shown  in  Fig.  B-l  .   Computed  values  are  listed  in 
Table  B-l. 

Once  the  function  %    is  known,  the  integral  I   may  be  rewritten  as 
follows . 


I2=i/^ 

Hir 


I  A!  J 


(-16) 


Similarly 


t-lR 


lj  JJrP}   1-7  /n2n,2 


* 


7r'  (JTfi'^dr'dz' 


(-17) 


,2 


The  scale  parameter  A   is  now  found  from  (-6) .   Notice  that  in 
computing  A   from  (-6)  the  factors  tt  /—  in  (-16)  and  (-17)  cancel  from  the 
result. 

Eqs.  (-16)  and  (-17)  illustrate  how  the  characteristic  integrals 
reduce  from  volume  to  surface  integrals  for  the  ax i- symmetric  case.   Eos. 
(-11)  and  (-13)  define.  ^  which  is  a  function  characteristic  of  axi-symmetric 
mean  flows. 
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O 

U. 
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UJ 
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O 
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TABLE  B-l 

COMPUTED  VALUES  OF  THE  ij<  FUNCTION 

(Refer  to  Eqs.  B-ll,  B-13) 

0.0  0.0000 

0.1  0.3475 

0.2  0.6555 

0.3  0.8955 

0.4  1.0562 

0.5  1.1433 

0.6  1.1739 

0.7  1.1687 

0.8  1.1461 

0.9  1.1187 

1.0  1.0936 

1.2  1.0581 

1.4  1.0387 
1.6  1.0279 
1.8  1.0213 
2.0  1.0169 

2.5  1.0105 
3.0  1.0072 
4.0  1.0040 
5.0  1.0025 

1.0000 
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APPENDIX  C 
Experimental  Data 

Various  useful  reference  items  of  experimental  data  are  summarized  in 
Tables  C-l  through  C-5  and  in  Figs.  C-l  and  C-2.   These  tabulations  and 
diagrams  are  largely  self-explanatory. 

Table  C-2  is  of  special  interest  in  that  it  shows  a  comparison 
between  theoretical  values  of  dissipation,  as  computed  from  the  heuristic 
model,  and  the  experimental  measurements  of  Laufer,  Ref.  (10).   The  agree- 
ment is  satisfactory. 

Laufer 's  measured  values  are  known  to  be  too  low  and  he  attempted  to 
extimate  better  values  as  shown. 

Shearing  rates  were  estimated  on  the  basis  of  a  theoretical  mean 
velocity  profile  using  Prandtl's  mixing  length  theory  and  Nikuradse's 
empirical  curve  for  the  mixing  length  I ,  namely, 

I   =  0.14  -  0.08r2  -  0.06r4  (C-l) 

The  information  in  Table  C-2  is  also  displayed  graphically  in  Figs. 
C-l  and  C-2.   There  is  no  particular  significance  to  the  fact  that  the 
ordinate  in  Fig  C-2  is  t**  rather  than  3  itself;  this  merely  reflects  the 
circumstance  that  a  slightly  different  notation  was  being  employed  at  the 
time  these  calculations  were  actually  made  than  that  which  is  used  in  this 
report. 

In  view  of  the  uncertainties  in  the  experimental  data  itself,  the 
agreement  between  theory  and  experiment  as  revealed  in  Figs.  C-l  and  C-2 
is  regarded  as  satisfactory. 


108 


The  remaining  tables  require  no  particular  comment  except,  perhaps, 
to  note  that  the  available  data  is  meagre  and  not  very  satisfactory.   Note 
that  the  data  given  is  mostly  one  dimensional,  with  slight  information  in 
the  second  dimension  and  none  at  all  in  the  third.   There  is  clearly  a 
need  for  more  experimental  work. 
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TABLE  C-l 

TURBULENT  ENERGY  DATA  FOR  A  PIPE 
Measurements  by  Laufer,  Ref.  (10) 

All  quantities  non-dimensionalized  using  radius  a  =  1  and  friction 
velocity  v*  =  1 . 


R  =  500,000 
e 


2E 


EH 


2E 


R  =  50,000 

e 


hi 


2lij  ij 


0 

0.1 

0.2 

0.3 

0.4 

0.5 

0.6 

0.7 

0.8 

0.9 

0.98 


1.65 

2.1 

1.80 

2.2 

2.05 

2.4 

2.45 

3.1 

2.90 

3.8 

3.50 

4.8 

A.  20 

6.2 

5.00 

8.7 

6.00 

13.0 

7.40 

22.8 

1.69 

2.15 

— 

1.80 

2.26 

— 

2.06 

2.72 

— 

2.41 

3.28 

1.5 

2.85 

3.96 

2.5 

3.48 

4.76 

3.7 

4.01 

6.23 

5.3 

4.49 

8.10 

7.6 

5.13 

11.3 

11.7 

5.66 

20.9 

22.0 

7.85(peak)   — 
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TABLE  C-2 

MEASURED  AND  COMPUTED  ENERGY  DISSIPATION 

RATES  IN  A  PIPE  AT  R  =  500,000 
__e * 

Measurements  by  Laufer,  Ref.  (10) 

Al  quantities  non-dimensionalized  using  radius  a  =  1  and  friction 
velocity  v*  =  1. 


(2E)7/6     (fifi')173 


-1 


Table  C-l     Calc.     Eq .  (1.3-18)    Eq .  (2.8-14)    Meas .     Est, 


0 

1.78 

2.94 

3.70 

1.43 

1.0 

2.0 

0.1 

2.00 

2.98 

3.70 

1.61 

0.2 

2.45 

3.09 

4.58 

1.65 

1.2 

2.5 

0.3 

2.95 

3.29 

6.71 

1.45 

0.4 

3.45 

3.62 

7.35 

1.71 

2.0 

4.0 

0.5 

4.30 

4.16 

7.40 

2.42 

0.6 

5.30 

5.03 

7.40 

3.61 

3.0 

6.0 

0.7 

6.70 

6.52 

7.40 

5.91 

0.8 

8.20 

9.58 

7.40 

10.60 

7.0 

13.5 

0.9 

10.00 

18.80 

7.40 

25.40 

15.0 

25.0 

Notes:   1.   Laufer  regards  his  measured  dissipation  data  as  too  low  and 
attempts  to  estimate  better  balues  as  indicated. 

1/3 
2.   The  quantity  (ftfi1)    was  obtained  from  the  theoretical 

velocity  profile  based  on  the  Prandtl-Nikuradse  mixing  length 

theory  using  Eq .  (C-l)  Refer  also  to  Eq .  (1.3-7). 
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TABLE  C-3 


AXIAL  VELOCITIES  IK  A  COLD  TURBULENT  JET 


AT  AXIAL  STATION  -|  =  20 


Measurement  by  Corrsin  and  Uberoi,  Ref.  (7) 


=  mean  axial  velocity 


V    =  axial  velocity  at  axis 
o 


v    =  root  mean  square  axial  perturbation 

V    1 

b,  ,„  =  radius  to  point  where  —  =  -= 


z    =  axial  coordinate 
D    =  diameter  of  jet 


'1/2 


'v  ^ 


0 

0.2 

0.4 

0.6 

0.8 

1.0 

1.2 

1.4 

1.6 

1.8 

2.0 

2.2 

2.4 

2.6 


1.0 

0.95 

0.85 

0.74 

0.625 

0.500 

0.40 

0.32 

0.25 

0.18 

0.135 

0.08 

0.055 

0.04 


0.205 

0.21 

0.215 

0.215 

0.20 

0.18 

0.165 

0.145 

0.12 

0.09 

0.07 

0.045 

0.03 

0.02 
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TABLE  C-4 


AXIAL  AND  RADIAL  VELOCITIES  IN  A  HEATED  TURBULENT  JET 


AT  AXIAL  STATION  ^  -  15 


Measurements  by  Corrsin  and  Uberoi,  Ref.  (7) 


u  =  root  mean  square  radial  perturbation 
Other  notation  as  in  Table  C-3 


1/2 


oJ 


T 


oJ 


oJ 


0 

0.2 

0.4 

0.6 

0.8 

1.0 

1.2 

1.4 

1.6 

1.8 

2.0 

2.2 

2.4 


1.0 

0.97 

0.88 

0.76 

0.62 

0.50 

0.38 

0.28 

0.205 

0.145 

0.10 

0.055 

0.025 


0.17 

0.17 

0.17 

0.17 

0.17 

0.15 

0.13 

0.10 

0.075 

0.05 

0.03 

0.02 

0.01 


0.175 

0.17 

0.16 

0.14 

0.12 

0.09 

0.075 

0.05 

0.035 

0.022 

0.013 

0.005 

0 
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TABLE  C-5 

AXIAL  VELOCITIES  IN  A  TURBULENT  JET 

AT  TWO  AXIAL  STATIONS 
Measurements  by  Laurence,  P.ef.  (8) 
Notation  as  in  Tables  C-3  and  C-A 


At 


=  16 


*  =  r 

r 

0 

0.5 

1.0 

1.5 

2.0 

2.5 

3.0 

3.5 

A.O 

A. 5 

5.0 

6.0 

7.0 

8.0 


0.A3 


0.A0 


0.31 


0.22 


0.15 


0.07 


0.082 
0.082 
0.082 
0.08A 
0.089 
0.091 
0.090 
0.08A 
0.07A 
0.061 
0.0A7 


At 

z 

D  = 

20 

fv  1 

*■  O' 

*■  O'' 

0.36 

0.073 

— 

0.073 

0.3A 

0.073 

— 

0.07A 

0.30 

0.075 

— 

0.075 

0.25 

0.075 

— 

0.075 

0.19 

0.075 

— 

O.C72 

o.i: 

35 

0.067 

0.08 

0.025 

0 
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